function setSelected2in3= phase1OnlySelect2(Y,X)
    
% v1.2 HB 6Mar2017 Now calls regress_hsb

[n,d]=size(X);
p=nthroot(d,3);
cube=reshape(1:d,p,p,p);
    cubeSelect=zeros(p,p,p); %initialisation for a cube which will count how many times each variable is selected as we traverse through the cube. Each element takes values {0,1,2,3};
    % traverse the cube
    for indL=1:p %layer index for traversal through cube
        for indR=1:p %row index for traversal through square in layer indL
            if size(nonzeros(cube(indR,:,indL)),1)>0
                subsetX=X(:,nonzeros(cube(indR,:,indL)));
                regressX=[ones(n,1),subsetX];
                    [~,tStat]=regress_hsb(Y,regressX);
                    tStat(1)=[];
                    [~,sortIndex]=sort(abs(tStat),'descend');
                    if length(sortIndex)>1
                        idxSelected=sortIndex(1:2); % only select the two most significant variables
                    else
                        idxSelected=sortIndex;
                    end
                    if size(idxSelected,1)>0
                        cubeSelect(indR,idxSelected,indL)=cubeSelect(indR,idxSelected,indL)+1;
                    end
            end
        end %indR
        for indC=1:p %column index for traversal through square in layer indL
            if size(nonzeros(cube(:,indC,indL)),1)>0
                subsetX=X(:,nonzeros(cube(:,indC,indL)));
                regressX=[ones(n,1),subsetX];
                [~,tStat]=regress_hsb(Y,regressX);
                    tStat(1)=[];
                    [~,sortIndex]=sort(abs(tStat),'descend');
                if length(sortIndex)>1
                    idxSelected=sortIndex(1:2); % only select the two most significant variables
                else
                    idxSelected=sortIndex;
                end
                if size(idxSelected,1)>0
                    cubeSelect(idxSelected,indC,indL)=cubeSelect(idxSelected,indC,indL)+1;
                end
            end
        end % indC
    end %indL
    % now traverse the cube in the layer direction
    for indR=1:p
        for indC=1:p
            if size(nonzeros(cube(indR,indC,:)),1)>0
                subsetX=X(:,nonzeros(cube(indR,indC,:)));
                regressX=[ones(n,1),subsetX];
                [~,tStat]=regress_hsb(Y,regressX);
                    tStat(1)=[];
                    [~,sortIndex]=sort(abs(tStat),'descend');
                if length(sortIndex)>1
                    idxSelected=sortIndex(1:2); % only select the two most significant variables
                else
                    idxSelected=sortIndex;
                end
                    if size(idxSelected,1)>0
                        cubeSelect(indR,indC,idxSelected)=cubeSelect(indR,indC,idxSelected)+1;
                    end
            end
        end % indC
    end %indR
    setSelected2in3=cube(find(cubeSelect>1)); % pick out the indices of the variables selected twice
