plskern

src/plskern.jl

@@ -103,8 +103,10 @@ end
 
 function plskern!(X::Matrix, Y::Union{Matrix, BitMatrix}, weights::Weight; kwargs...)
     par = recovkw(ParPlsr, kwargs).par
+    ## Specific for Plsda functions
     Q = eltype(X)
     isa(Y, BitMatrix) ? Y = convert.(Q, Y) : nothing
+    ## End
     n, p = size(X)
     q = nco(Y)
     nlv = min(n, p, maximum(par.nlv)) # 'maximum' required for plsravg 
@@ -121,9 +123,12 @@ function plskern!(X::Matrix, Y::Union{Matrix, BitMatrix}, weights::Weight; kwarg
         fcenter!(X, xmeans)
         fcenter!(Y, ymeans)
     end
-    D = Diagonal(weights.w)
-    XtY = X' * (D * Y)                   # = Xd' * Y = X' * D * Y  (Xd = D * X   Very costly!!)
-    #XtY = X' * (weights.w .* Y)         # Can create OutOfMemory errors for very large matrices
+    ## XtY 
+    fweight!(Y, weights.w)
+    XtY = X' * Y
+    ## Old
+    ## D = Diagonal(weights.w)
+    ## XtY = X' * (D * Y)    # Xd = D * X   Very costly!!
     ## Pre-allocation
     T = similar(X, n, nlv)
     W = similar(X, p, nlv)

src/plsnipals.jl

@@ -88,7 +88,7 @@ function plsnipals!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
         zp ./= tt
         mul!(c, Y', dt)
         c ./= tt                      
-        # deflation with respect to t (asymetric PLS)
+        # deflation with respect to t: asymetric PLS
         X .-= t * zp'
         Y .-= t * c'
         # end

src/plsrosa.jl

@@ -70,7 +70,7 @@ function plsrosa!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
     c = similar(X, q)
     # End
     @inbounds for a = 1:nlv
-        XtY .= X' * (D * Y)
+        XtY .= X' * D * Y
         if q == 1
             w .= vec(XtY)
             w ./= normv(w)
@@ -83,7 +83,7 @@ function plsrosa!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
             t .= t .- z * inv(z' * D * z) * z' * (D * t)
             z = vcol(W, 1:(a - 1))
             w = w .- z * (z' * w)
-            w ./= sqrt(dot(w, w))
+            w ./= normv(w)
         end
         dt .= weights.w .* t
         tt = dot(t, dt)

src/plsrout.jl

@@ -85,5 +85,5 @@ function plsrout!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
     w .*= wtal(d; a = quantile(d, 1 - par.prm))
     w .*= weights.w
     w[isequal.(w, 0)] .= 1e-10
-    plskern(X, Y, mweight(w); kwargs...)
+    plskern!(X, Y, mweight(w); kwargs...)
 end

src/plssimp.jl

@@ -53,8 +53,9 @@ function plssimp!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
         fcenter!(X, xmeans)
         fcenter!(Y, ymeans)
     end
-    D = Diagonal(weights.w)
-    XtY = X' * (D * Y)   
+    ## XtY 
+    fweight!(Y, weights.w)
+    XtY = X' * Y
     ## Pre-allocation
     T = similar(X, n, nlv)
     W = similar(X, p, nlv)

src/plswold.jl

@@ -45,7 +45,6 @@ function plswold!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
     n, p = size(X)
     q = nco(Y)
     nlv = min(par.nlv, n, p)
-    sqrtw = sqrt.(weights.w)
     xmeans = colmean(X, weights) 
     ymeans = colmean(Y, weights)   
     xscales = ones(Q, p)
@@ -60,8 +59,9 @@ function plswold!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
         fcenter!(Y, ymeans)
     end
     # Row metric
-    X .= sqrtw .* X
-    Y .= sqrtw .* Y
+    sqrtw = sqrt.(weights.w)
+    fweight!(X, sqrtw)
+    fweight!(Y, sqrtw)
     ## Pre-allocation
     Tx = similar(X, n, nlv)
     Wx = similar(X, p, nlv)
@@ -76,9 +76,9 @@ function plswold!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
     px   = copy(wx)
     niter = zeros(nlv)
     # End
-    @inbounds for a = 1:nlv
-        tx .= X[:, 1]
-        ty .= Y[:, 1]
+    @inbounds for a = 1:nlv       
+        tx .= vcol(X, 1)
+        ty .= vcol(Y, 1)
         cont = true
         iter = 1
         wx .= rand(p)
@@ -110,7 +110,7 @@ function plswold!(X::Matrix, Y::Matrix, weights::Weight; kwargs...)
         Wytild[:, a] .= wytild
         TTx[a] = ttx
     end
-    Tx .= (1 ./ sqrtw) .* Tx
+    fweight!(Tx, 1 ./ sqrtw) 
     Rx = Wx * inv(Px' * Wx)
     Plsr(Tx, Px, Rx, Wx, Wytild, TTx, xmeans, xscales, ymeans, yscales, weights, niter, par)
 end