diff --git a/src/PartialLeastSquaresRegressor.jl b/src/PartialLeastSquaresRegressor.jl index fa1591a..9cf19e9 100644 --- a/src/PartialLeastSquaresRegressor.jl +++ b/src/PartialLeastSquaresRegressor.jl @@ -9,6 +9,7 @@ include("types.jl") include("pls1.jl") include("pls2.jl") include("kpls.jl") +include("opls.jl") include("method.jl") include("mlj_interface.jl") diff --git a/src/opls.jl b/src/opls.jl new file mode 100644 index 0000000..59e62d6 --- /dev/null +++ b/src/opls.jl @@ -0,0 +1,53 @@ +# based on the code from https://github.com/BiRG/pyopls +# reference: Johan Trygg and Svante Wold. Orthogonal projections to latent structures (O-PLS). +# J. Chemometrics 2002; 16: 119-128. DOI: 10.1002/cem.695 + + + +##filtering and learning - Ortogonal PLS, one feature +function fitting(model::OPLS1Model, + X::AbstractArray{T}, + Y::Vector{T}) where T<:AbstractFloat + + w=X'*Y # calculate weight vector + w=w/norm(w) # normalization + + for i in 1:model.n_ortho_components + t=X*w # calculate scores vector nrows + p=X'*t/(t'*t) # calculate loadings of X ncols + wosc=p-(w'*p)/(w'*w)*w # orthogonal weight ncols + wosc=wosc/norm(wosc) # normalization ncols + tosc=X*wosc # orthogonal components nrows + posc=X'*tosc/(tosc'*tosc) # loadings ncols + + X=X-tosc*posc' # remove orthogonal components + + model.W_ortho[:,i]=wosc # weights orthogonal to y + model.P_ortho[:,i]=posc # loadings orthogonal to y + model.T_ortho[:,i]=tosc # scores orthogonal to y + end + + # X is now with orthogonal signal components removed + return X, model.T_ortho +end + +# remove orthogonal components from Xt +# return filtered Xt and weightings of orthogonal components +function filter!(model::OPLS1Model, + X::AbstractArray{T}) where T<:AbstractFloat + nrow = size(X,1) + # ortogonal weights + ortho = zeros(T, nrow, model.n_ortho_components) + + for i in 1:model.n_ortho_components + R = X * model.W_ortho[:,i] + X = X - R * model.P_ortho[:,i]' + ortho[:,i] = R + end + + return X, ortho +end + +function component(model::OPLS1Model{T},i) where T<:AbstractFloat + return model.W_ortho[:,i] +end \ No newline at end of file diff --git a/src/pls1.jl b/src/pls1.jl index d92f1dc..a08f0dc 100644 --- a/src/pls1.jl +++ b/src/pls1.jl @@ -40,3 +40,23 @@ function predictor(model::PLS1Model{T}, return Y end + +function filter!(model::PLS1Model{T}, + X::AbstractArray{T}) where T<:AbstractFloat + + nrows = size(X,1) + X_proj = zeros(T, nrows, model.nfactors) + Y_pred = zeros(size(X,1)) + + for i = 1:model.nfactors + X_proj[:,i] = X*model.W[:,i] + X = X - X_proj[:,i]*model.P[:,i]' + Y_pred = Y_pred + X_proj[:,i]*model.b[i] # building prediction + end + + return X, X_proj, Y_pred +end + +function component(model::PLS1Model{T},i) where T<:AbstractFloat + return model.W[:,i] +end \ No newline at end of file diff --git a/src/pls2.jl b/src/pls2.jl index df9be5f..0609178 100644 --- a/src/pls2.jl +++ b/src/pls2.jl @@ -61,3 +61,27 @@ function predictor(model::PLS2Model{T}, return Y end +function filter!(model::PLS2Model{T}, + X::AbstractArray{T}) where T<:AbstractFloat + + nrows = size(X,1) + X_proj = zeros(T, nrows, model.nfactors) + + W,Q,P = model.W,model.Q,model.P + nfactors = model.nfactors + Y = zeros(T,size(X,1),model.ntargetcols) + #println("nfactors: ",nfactors) + for i = 1:nfactors + #R = X*W[:,i] + X_proj[:,i] = X*model.W[:,i] + X = X - X_proj[:,i] * P[:,i]' + Y = Y + X_proj[:,i] * Q[:,i]' + end + + return X,X_proj,Y +end + + +function component(model::PLS2Model{T},i) where T<:AbstractFloat + return model.W[:,i] +end \ No newline at end of file diff --git a/src/types.jl b/src/types.jl index 60d1909..f525013 100644 --- a/src/types.jl +++ b/src/types.jl @@ -93,3 +93,30 @@ function PLSModel(X::Matrix{T}, kernel, width) end + + +################################################################################ +#### OPLS1 type +mutable struct OPLS1Model{T<:AbstractFloat} <:PLSModel{T} + W_ortho::Matrix{T} # : weights orthogonal to y + P_ortho::Matrix{T} # : loadings orthogonal to y + T_ortho::Matrix{T} # : scores orthogonal to y + + n_ortho_components::Int # +end + +## OPLS1: constructor +function OPLS1Model(X::Matrix{T}, + Y::Vector{T}, + n_components::Int) where T<:AbstractFloat + + (nrows,ncols) = size(X) + + ## Allocation + return OPLS1Model( + zeros(T, ncols, n_components), ## W_ortho + zeros(T, ncols, n_components), ## P_ortho + zeros(T, nrows, n_components), ## T_ortho + n_components ## n_components + ) +end \ No newline at end of file diff --git a/test/Project.toml b/test/Project.toml new file mode 100644 index 0000000..4524baa --- /dev/null +++ b/test/Project.toml @@ -0,0 +1,8 @@ +[deps] +MLJBase = "a7f614a8-145f-11e9-1d2a-a57a1082229d" +MLJModelInterface = "e80e1ace-859a-464e-9ed9-23947d8ae3ea" +Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +ScientificTypes = "321657f4-b219-11e9-178b-2701a2544e81" +ScientificTypesBase = "30f210dd-8aff-4c5f-94ba-8e64358c1161" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"