Add files via upload

change from gradually remove the Orthogonal variation to remove them together

Author: guokai8

R/o2pls.R

@@ -20,8 +20,10 @@
 #'    \item{CXosc}{Orthogonal \eqn{Y} weights}
 #'    \item{BU}{Regression coefficient in \code{Tt} ~ \code{U}}
 #'    \item{BT}{Regression coefficient in \code{U} ~ \code{Tt}}
-#'    \item{R2Xhat}{Prediction of \eqn{X} with \eqn{Y}}
-#'    \item{R2Yhat}{Prediction of \eqn{Y} with \eqn{X}}
+#'    \item{Xhat}{Prediction of \eqn{X} with \eqn{Y}}
+#'    \item{Yhat}{Prediction of \eqn{Y} with \eqn{X}}
+#'    \item{R2Xhat}{Variation of the predicted \eqn{X} as proportion of variation in \eqn{X}}
+#'    \item{R2Yhat}{Variation of the predicted \eqn{Y} as proportion of variation in \eqn{Y}}
 #'    \item{R2X}{Variation of the modeled part in \eqn{X} (defined by Joint + Orthogonal variation) as proportion of total variation in \eqn{X}}
 #'    \item{R2Y}{Variation of the modeled part in \eqn{Y} (defined by Joint + Orthogonal variation) as proportion of total variation in \eqn{Y}}
 #'    \item{R2Xcorr}{Variation of the joint part in \eqn{X} }
@@ -94,7 +96,6 @@ o2pls <- function(X,Y,nc,nx,ny,scale=FALSE,center=FALSE){
     CXosc <- matrix(0, ncy, nny)
 
     # Calculate the PCA components of Y'T
-   ## nn<-nc+max(nx+ny)
     ###borrow idea from OmicsPLS
     ns <- nc + max(nx, ny)
     CDW <- eigenmult(t(Y),X)
@@ -106,50 +107,43 @@ o2pls <- function(X,Y,nc,nx,ny,scale=FALSE,center=FALSE){
 #    Yscore <- Y %*% Yloading
     Yscore <- eigenmult(Y,Yloading)
     if(nx>0){
-        for(i in 1:nx){
+        #####Notice the change #### Now remove the orthogonal by doing once
+        #### Still thinking to do it one by one #######
+        #### In previous version I did this one by one #######
           #  Exy <- X - Xscore%*%t(Xloading)
             Exy <- X - eigenmult(Xscore, t(Xloading))
           #  wosc<- svd(t(Exy)%*%Xscore,nu=1,nv=0)$u
-            wosc<- svd(eigenmult(t(Exy),Xscore),nu=1,nv=0)$u
+            wosc<- svd(eigenmult(t(Exy),Xscore),nu=nx,nv=0)$u
            # t_yo <- X%*%wosc
             t_yo <- eigenmult(X,wosc)
          #   p_yo <- t(solve(t(t_yo)%*%t_yo)%*%t(t_yo)%*%X)
             p_yo <- t(eigenthree(solve(AtA(t_yo)),t(t_yo),X))
           #  X <- X - t_yo%*%t(p_yo)
             X <- X - eigenmult(t_yo,t(p_yo))
-            TYosc[,i] <- t_yo
-            PYosc[,i] <- p_yo
-            WYosc[,i] <- wosc	
-            #Xscore <- X%*%Xloading
-            CDW <- eigenmult(t(Y), X)
-            cdw <- svd(CDW, nu = nc, nv = nc)
-            Xloading <- cdw$v
+            TYosc <- t_yo
+            PYosc <- p_yo
+            WYosc <- wosc	
             Xscore <- eigenmult(X, Xloading)
-        }
+ 
     }
     if(ny>0){
-        for(i in 1:ny){
            # Fxy <- Y - Yscore %*% t(Yloading)
             Fxy <- Y - eigenmult(Yscore, t(Yloading))
            # cosc <- svd(t(Fxy)%*%Yscore,nu=1,nv=0)$u
-            cosc <- svd(eigenmult(t(Fxy), Yscore),nu=1,nv=0)$u
+            cosc <- svd(eigenmult(t(Fxy), Yscore),nu=ny,nv=0)$u
            # u_xo <- Y%*%cosc
             u_xo <- eigenmult(Y, cosc)
            # p_xo <- t(solve(t(u_xo)%*%u_xo)%*%t(u_xo)%*%Y)
             p_xo <- t(eigenthree(solve(AtA(u_xo)),t(u_xo),Y))
           #  Y <- Y - u_xo%*%t(p_xo)
             Y <- Y - eigenmult(u_xo, t(p_xo))
-            UXosc[,i] <- u_xo
-            PXosc[,i] <- p_xo
-            CXosc[,i] <- cosc	
+            UXosc <- u_xo
+            PXosc <- p_xo
+            CXosc <- cosc	
           #  Yscore <- Y%*%Yloading
-            CDW <- eigenmult(t(Y), X)
-            cdw <- svd(CDW, nu = nc, nv = nc)
-            Yloading <- cdw$u
             Yscore <- eigenmult(Y, Yloading)
-        }
     }
-    ##calculate again
+    ##calculate the score and loading again
     CDW <- eigenmult(t(Y), X)
     cdw <- svd(CDW, nu = nc, nv = nc)
     Xloading <- cdw$v
@@ -168,8 +162,13 @@ o2pls <- function(X,Y,nc,nx,ny,scale=FALSE,center=FALSE){
     R2Y <- 1- s2(Fxy)/SSY
     Xh <- eigenthree(Yscore, BU, t(Xloading))
     Yh <- eigenthree(Xscore, BT, t(Yloading))
-    R2Xhat <- 1 - s2(Xh-Xt)/SSX
-    R2Yhat <- 1 - s2(Yh-Yt)/SSY
+    #####
+   # R2Xhat <- 1 - s2(Xh-Xt)/SSX
+   # R2Yhat <- 1 - s2(Yh-Yt)/SSY
+    #### #########
+    R2Xhat <- (s2(Yscore %*% BU)/s2(Xt))
+    R2Yhat <- (s2(Xscore %*% BT)/s2(Yt))
+    ####
     R2Xp <- (s2(Xh)/SSX)/R2Xcorr
     R2Yp <- (s2(Yh)/SSY)/R2Ycorr
     Qx <- Q(Xt,Xh)
@@ -178,10 +177,10 @@ o2pls <- function(X,Y,nc,nx,ny,scale=FALSE,center=FALSE){
     varYj = apply(Yscore,2,function(x)sum(x^2))/SSY
     varXorth = apply(PYosc,2,function(x)sum(x^2))*apply(TYosc,2,function(x)sum(x^2))/SSX
     varYorth = apply(PXosc,2,function(x)sum(x^2))*apply(UXosc,2,function(x)sum(x^2))/SSY
-    rownames(Xscore) <- rownames(TYosc) <- rownames(Exy) <- rownames(Xt)
-    rownames(Yscore) <- rownames(UXosc) <- rownames(Fxy) <- rownames(Yt)
-    rownames(Xloading) <- rownames(PYosc) <- rownames(WYosc) <- colnames(Exy) <- colnames(Xt)
-    rownames(Yloading) <- rownames(PXosc) <- rownames(CXosc) <- colnames(Fxy) <- colnames(Yt)
+    rownames(Xscore) <- rownames(TYosc) <- rownames(Exy) <- rownames(Xh) <- rownames(Xt)
+    rownames(Yscore) <- rownames(UXosc) <- rownames(Fxy) <- rownames(Yh) <- rownames(Yt)
+    rownames(Xloading) <- rownames(PYosc) <- rownames(WYosc) <- colnames(Exy) <- colnames(Xh)<- colnames(Xt)
+    rownames(Yloading) <- rownames(PXosc) <- rownames(CXosc) <- colnames(Fxy) <- colnames(Yh)<-colnames(Yt)
     colnames(Xscore) <- colnames(Yscore) <- colnames(Xloading) <- colnames(Yloading) <- paste0("LV",1:nc)
     colnames(TYosc) <- colnames(PYosc) <- colnames(WYosc) <- paste0("LV", 1:nnx)
     colnames(UXosc) <- colnames(PXosc) <- colnames(CXosc) <- paste0("LV", 1:nny)
@@ -194,6 +193,7 @@ o2pls <- function(X,Y,nc,nx,ny,scale=FALSE,center=FALSE){
                 PXosc = PXosc, PYosc = PYosc,
                 WYosc = WYosc, CXosc = CXosc,
                 R2Xhat = R2Xhat, R2Yhat = R2Yhat,
+                Xhat = Xh, Yhat= Yh,
                 Qx = Qx, Qy = Qy,
                 varXj = varXj, varYj = varYj, 
                 varXorth = varXorth, varYorth =varYorth,