This commit is contained in:
Joakim Skogholt 2023-05-10 22:14:23 +02:00
parent 239b021eef
commit 24b0117d48

View file

@ -1,35 +1,46 @@
""" """
function PCR(X, y, kmax, centre=True, standardize=true) function PCR(X, y, kmax, centre=true)#, standardize=true)
Principal Component Regression (PCR). Principal Component Regression (PCR).
Inputs: Data matrix, response vector, maximum number of components. Inputs: Data matrix, response vector, maximum number of components.
A constant term is included in the modeling. A constant term is included in the modeling if centre==true.
Outputs: B (matrix of size (p+1) x kmax), U, s, V Outputs: B (matrix of size (p+1) x kmax), U, s, V
ADD STANDARDIZATION?? (NEED TO THINK THROUGH PREDICTION WITH NEW DATA)
X, y = importData("Beer"); X, y = importData("Beer");
B, \\_ = PCR(X, y, 10, true, false); B, \\_ = PCR(X, y, 10, true, false);
""" """
function PCR(X, y, kmax, centre=true, standardize=true) function PCR(X, y, kmax, centre::Bool=true)#, standardize=true)
stdX = std(X, dims=1);
mX = mean(X, dims=1); mX = mean(X, dims=1);
my = mean(y, dims=1); my = mean(y, dims=1);
stdX = std(X, dims=1); y = y .- my;
if centre if centre
X = X .- mX; X = X .- mX;
end end
if standardize #if standardize
X = X ./ stdX; # X = X ./ stdX;
end #end
U, s, V = svd(X, full=false); U, s, V = svd(X, full=false);
q = s[1:kmax].^(-1) .*(U[:,1:kmax]'y); q = s[1:kmax].^(-1) .*(U[:,1:kmax]'y);
B = cumsum(V[:,1:kmax] .* q', dims=2); B = cumsum(V[:,1:kmax] .* q', dims=2);
b0 = my .- mX * B
B = [b0; B]; if centre
b0 = my .- mX * B
B = [b0; B];
end
return B, U, s, V return B, U, s, V
end end