modldemo.m

偏最小二乘算法在MATLAB中的实现

echo on
%MODLDEMO Demonstrates MODLMKER and MODLPRED functions
%
% This demonstration illustrates the use of the MODLMKER
% function in the PLS_Toolbox for general regression
% modelling. The MODLRDER function for reading models and
% MODLPRED function for making new predictions are also
% demonstrated.

echo off
%  Copyright
%  Barry M. Wise
%  July 1994
echo on

% The data we are going to work with (once again) is from a 
% Liquid-Fed Ceramic Melter (LFCM).  We will develop a 
% correlation between temperatures in the molten glass tank 
% and the tank level.

% Lets start by loading and plotting the data.  Hit a key when
% you are ready.
pause

echo off
load plsdata
subplot(2,1,1)
plot(xblock1);
title('X-block Data (Predictor Variables) for PLS Demo');
xlabel('Sample Number');
ylabel('Temperature (C)');
subplot(2,1,2)
plot(yblock1)
title('Y-Block Data (Predicted Variable) for PLS Demo');
xlabel('Sample Number');
ylabel('Level (Inches)');
echo on

% You can probably already see that there is a very regular
% variation in the temperature data and that it appears to
% correlate with the level data. This is because there is
% steep temperature gradient in the molten glass, and when
% the level changes, glasses of different temperatures pass
% by the location of the thermocouples.
pause

% Lets use the fact that temperature correlates with
% level to build a PLS, PCR or RR model that uses temperature
% to predict level. We will use the MODLMKER function, which
% will prompt you for your choices of data scaling, regression
% technique and cross validation method. For this data set,
% mean centering and cross validation using contiguous blocks
% of data is most appropriate. You may want to run this several
% times and see how it turns out when different options are
% chosen.

% We will now call the MODLMKER function. Note that the only
% inputs are the predictor and predicted variable matrices.
echo off
subplot(1,1,1)
echo on

pause

model = modlmker(xblock1,yblock1);

pause

% Now that you have developed a model, we can use MODLRDER
% to "read" it. 

modlrder(model)

% Note how the model you created is stored with a date and
% time stamp and some information about how the model was
% created. The final regression vector and scaling parameters
% are also saved with the final model.

pause

% The model created by MODLMKER can now be used with a new
% data set to make new predictions as follows. 

ypred = modlpred(xblock2,model);

% We can plot up the actual and predicted values and see
% how we did.

pause
echo off
plot(1:200,yblock2,'-y',1:200,yblock2,'og'), hold on
plot(1:200,ypred,'-r',1:200,ypred,'+r')
xlabel('Sample Number (time)')
ylabel('Tank Level (inches)')
title('Actual (o) and Predicted (+) Tank Level')

pause
echo on
% From earlier experience, we know that the best PLS model
% of this system has 5 LVs. Lets quickly build that PLS model
% and see how it compares to the model you just developed
% with MODLMKER. We will also calculate the MLR model for
% comparison.

[mxblock1,mx] = mncn(xblock1);
[myblock1,my] = mncn(yblock1);
[p,q,w,t,u,b,ssqdif] = pls(mxblock1,myblock1,5);
mlrmod = mxblock1\myblock1;
sxblock2 = scale(xblock2,mx);
ypredopt = plspred(sxblock2,b,p,q,w,5);
ypredmlr = sxblock2*mlrmod;
sypredopt = rescale(ypredopt,my);
sypredmlr = rescale(ypredmlr,my);
pause

% We can now add the PLS prediction to the figure.
echo off
plot(1:200,sypredopt,'-b',1:200,sypredopt,'*b'), hold off
echo on

% Hopefully, the predictions that you came up with using
% MODLMKER and MODLPRED aren't too much different from the
% PLS predictions we just added to the plot.
pause

% We can calculate the root mean square error of prediction 
% for your model and the PLS model as follows:

echo off
plsssq = sqrt(sum((yblock2-sypredopt).^2)/200);
modlssq = sqrt(sum((yblock2-ypred).^2)/200);
mlrssq = sqrt(sum((yblock2-sypredmlr).^2)/200);

disp('  PLS error Your Model MLR error'), 
disp([plsssq modlssq mlrssq])
echo on

% So here we see that the PLS model is slightly
% better than the MLR model. How is your model?