R/d.spls.GL.R

In dual.spls: Dual Sparse Partial Least Squares Regression

#' Dual Sparse Partial Least Squares (Dual-SPLS) regression for the group lasso norms
#' @description
#' The function performs dimensional reduction with the group lasso norms. Three norms are available
#' where \code{G} is the number of groups, the vectors \eqn{w_g} hold the coordinates of \eqn{w}
#' for the observations belonging to the group \eqn{g} and \eqn{\alpha_g}, \eqn{\lambda_g} and \eqn{\gamma_g} are all positive scalars.
#' \itemize{
#' \item Norm A *(generalized norm)*: \eqn{\Omega_g(w)=\|w_g\|_2+ \lambda_g \|w_g\|_1} where
#' \eqn{\Omega(w)=\sum_{g} \alpha_g \Omega_g(w)=1 \textrm{ and } \sum_{g=1}^G \alpha_g=1},
#
#' \item Norm B *(particular case)*: \eqn{\Omega(w)=\|w\|_2+\sum_{g=1}^G \lambda_g\|w_g\|_1},
#'
#' \item Norm C *(particular case)*: \eqn{\Omega(w)=\sum_{g=1}^G \alpha_g \|w \|_2+\sum_{g=1}^G \lambda_g \|w_g \|_1}
#'  where
#' \eqn{\sum_{g=1}^G \alpha_g=\sum_{g=1}^G \gamma_g=1} \cr and \eqn{\Omega(w_g)=\gamma_g}.
#' }
#'
#' Dual-SPLS for the group lasso norms has been designed to confront the situations where the predictors
#' variables can be divided into distinct meaningful groups. Each group is constrained by an independent
#' threshold as in the dual sparse lasso methodology,
#' that is each \eqn{w_g} will be collinear to a vector \eqn{z_{\nu_g}} built from the coordinate of \eqn{z}
#' and constrained by the threshold \eqn{\nu_g}.
#'
#' Three variants are defined here depending on the groups combination in the global norm and the weights
#' assigned to each group. They all give the same result as the lasso norm for \eqn{G=1},
#' \itemize{
#' \item Norm A is the generalized norm of the group lasso. applies the lasso norm for each group individually while constraining the overall norm. Moreover,
#' the Euclidean norm of each \eqn{w_g} is computed while minimizing the root mean squares error of prediction,
#' \item Norm B is a particular case and a genuine alternative similar to the lasso-like norm,
#' \item Norm C is another particular case that assigns user to define weights for each group.
#' }
#' @usage d.spls.GL(X,y,ncp,ppnu,indG,gamma=NULL,norm="A",verbose=FALSE)
#' @param X a numeric matrix of predictors values of dimension \code{(n,p)}. Each row represents one observation and each column one predictor variable.
#' @param y a numeric vector or a one column matrix of responses. It represents the response variable for each observation.
#' @param ncp a positive integer. \code{ncp} is the number of Dual-SPLS components.
#' @param ppnu a positive real value or a vector of length the number of groups, in \eqn{[0,1]}.
#' \code{ppnu} is the desired proportion of variables to shrink to zero for each component and for each group.
#' @param indG a numeric vector of group index for each observation.
#' @param gamma a numeric vector of the norm \eqn{\Omega} of each \eqn{w_g} in case \code{norm="C"}.
#' @param norm a character specifying the norm chosen between A, B and C. Default value is \code{A}.
#' @param verbose a Boolean value indicating whether or not to display the iterations steps. Default value is \code{FALSE}.
#' @details
#' The resulting solution for \eqn{w} and hence for the coefficients vector, in the case of \code{d.spls.GL}, has
#' a simple closed form expression (ref) deriving from the fact that for each group \eqn{g}, \eqn{w_g}
#' is collinear to a vector
#' \deqn{z_{\nu,g}=\textrm{sign}({z_g})(|z_g|-\nu_g)_+.}
#' Here, for each group \eqn{g}, \eqn{\nu_g} is the threshold for which \code{ppnu} of the group \eqn{g} of
#' the absolute values of the coordinates of \eqn{z_j} are greater than \eqn{\nu_g}. The norms differ in the value of the threshold for each group,
#' that is the expression of \eqn{\nu_g}. (see reference for detail)
#' @return A \code{list} of the following attributes
#' \item{Xmean}{the mean vector of the predictors matrix \code{X}.}
#' \item{scores}{the matrix of dimension \code{(n,ncp)} where \code{n} is the number of observations. The \code{scores} represents
#' the observations in the new component basis computed by the compression step
#' of the Dual-SPLS.}
#' \item{loadings}{the matrix of dimension \code{(p,ncp)} that represents the Dual-SPLS components.}
#' \item{Bhat}{the matrix of dimension \code{(p,ncp)} that regroups the regression coefficients for each component.}
#' \item{intercept}{the vector of length \code{ncp} representing the intercept values for each component.}
#' \item{fitted.values}{the matrix of dimension \code{(n,ncp)} that represents the predicted values of \code{y}}
#' \item{residuals}{the matrix of dimension \code{(n,ncp)} that represents the residuals corresponding
#' to the difference between the responses and the fitted values.}
#' \item{lambda}{the matrix of dimension \code{(G,ncp)} collecting the parameters of sparsity \eqn{\lambda_g} used to fit the model at each iteration and for each group.}
#' \item{alpha}{the matrix of dimension \code{(G,ncp)} collecting the constraint parameters \eqn{\alpha_g}  used to fit the model at each iteration and for each group when the norm chosen is \code{B} or \code{C}.}
#' \item{zerovar}{the matrix of dimension \code{(G,ncp)} representing the number of variables shrank to zero per component and per group.}
#' \item{PP}{the vector of length \code{G} specifying the number of variables in each group.}
#' \item{ind_diff0}{the list of \code{ncp} elements representing the index of the none null regression coefficients elements.}
#' \item{type}{a character specifying the Dual-SPLS norm used. In this case it is either \code{GLA}, \code{GLB} or \code{GLC}. }
#' @author Louna Alsouki François Wahl
#' @seealso [dual.spls::d.spls.GLA], [dual.spls::d.spls.GLB], [dual.spls::d.spls.GLC]
#'
#'
#' @examples
#' ### load dual.spls library
#' library(dual.spls)
#' oldpar <- par(no.readonly = TRUE)
#'
#' ####two predictors matrix
#' ### parameters
#' n <- 100
#' p <- c(50,100)
#' nondes <- c(20,30)
#' sigmaondes <- c(0.05,0.02)
#' data=d.spls.simulate(n=n,p=p,nondes=nondes,sigmaondes=sigmaondes)
#'
#' X <- data$X
#' X1 <- X[,(1:p[1])]
#' X2 <- X[,(p[1]+1):sum(p)]
#' y <- data$y
#'
#' indG <-c(rep(1,p[1]),rep(2,p[2]))
#'
#' #fitting the model
#' ncp <- 10
#' ppnu <- c(0.99,0.9)
#'
#' # norm A
#' mod.dsplsA <- d.spls.GL(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,norm="A",verbose=TRUE)
#' n <- dim(X)[1]
#' p <- dim(X)[2]
#'
#' str(mod.dsplsA)
#'
#' ### plotting the observed values VS predicted values
#' plot(y,mod.dsplsA$fitted.values[,6], xlab="Observed values", ylab="Predicted values",
#'  main="Observed VS Predicted for 6 components")
#' points(-1000:1000,-1000:1000,type='l')
#'
#' ### plotting the regression coefficients
#'
#' i=6
#' nz=mod.dsplsA$zerovar[,i]
#' plot(1:dim(X)[2],mod.dsplsA$Bhat[,i],type='l',
#'     main=paste(" Dual-SPLS (GLA), ncp =", i, " #0coef =", nz[1], "/", dim(X1)[2]
#'     , " #0coef =", nz[2], "/", dim(X2)[2]),
#'     ylab='',xlab='' )
#' inonz=which(mod.dsplsA$Bhat[,i]!=0)
#' points(inonz,mod.dsplsA$Bhat[inonz,i],col='red',pch=19,cex=0.5)
#' legend("topright", legend ="non null values", bty = "n", cex = 0.8, col = "red",pch=19)
#'
#' # norm B
#' mod.dsplsB <- d.spls.GL(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,norm="B",verbose=TRUE)
#'
#' str(mod.dsplsB)
#'
#' ### plotting the observed values VS predicted values
#' plot(y,mod.dsplsB$fitted.values[,6], xlab="Observed values", ylab="Predicted values",
#' main="Observed VS Predicted for 6 components")
#' points(-1000:1000,-1000:1000,type='l')
#'
#' ### plotting the regression coefficients
#'
#' i=6
#' nz=mod.dsplsB$zerovar[,i]
#' plot(1:dim(X)[2],mod.dsplsB$Bhat[,i],type='l',
#'     main=paste(" Dual-SPLS (GLB), ncp =", i, " #0coef =", nz[1], "/", dim(X1)[2]
#'     , " #0coef =", nz[2], "/", dim(X2)[2]),
#'     ylab='',xlab='' )
#' inonz=which(mod.dsplsB$Bhat[,i]!=0)
#' points(inonz,mod.dsplsB$Bhat[inonz,i],col='red',pch=19,cex=0.5)
#'
#' legend("topright", legend ="non null values", bty = "n", cex = 0.8, col = "red",pch=19)
#'
#' # norm C
#' mod.dsplsC <- d.spls.GL(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,gamma=c(0.5,0.5),norm="C",verbose=TRUE)
#' n <- dim(X)[1]
#' p <- dim(X)[2]
#'
#' str(mod.dsplsC)
#'
#' ### plotting the observed values VS predicted values
#' plot(y,mod.dsplsC$fitted.values[,6], xlab="Observed values", ylab="Predicted values",
#' main="Observed VS Predicted for 6 components")
#' points(-1000:1000,-1000:1000,type='l')
#'
#' ### plotting the regression coefficients
#'
#' i=6
#' nz=mod.dsplsC$zerovar[,i]
#' plot(1:dim(X)[2],mod.dsplsC$Bhat[,i],type='l',
#'     main=paste(" Dual-SPLS (GLC), ncp =", i, " #0coef =", nz[1], "/", dim(X1)[2]
#'     , " #0coef =", nz[2], "/", dim(X2)[2]),
#'     ylab='',xlab='' )
#' inonz=which(mod.dsplsC$Bhat[,i]!=0)
#' points(inonz,mod.dsplsC$Bhat[inonz,i],col='red',pch=19,cex=0.5)
#' legend("topright", legend ="non null values", bty = "n", cex = 0.8, col = "red",pch=19)
#'
#' par(oldpar)
#' @export


d.spls.GL<- function(X,y,ncp,ppnu,indG,gamma=NULL,norm="A",verbose=FALSE)
{
  if (norm=="A")
  {
    mod.dspls=d.spls.GLA(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,verbose=verbose)
  }
  if (norm=="B")
  {
    mod.dspls=d.spls.GLB(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,verbose=verbose)
  }

  if (norm=="C")
  {
    mod.dspls=d.spls.GLC(X=X,y=y,ncp=ncp,ppnu=ppnu,indG=indG,gamma=gamma,verbose=verbose)
  }
  return(mod.dspls)
}