www.gusucode.com > mbctools 工具箱 matlab 源码程序 > mbctools/@mdev_local/twostage.m

    function TSModels= twostage(mdev,selrf)
%TWOSTAGE make twostage models for mdev_local
%
% TSModels= twostage(mdev,selrf)

%  Copyright 2000-2015 The MathWorks, Inc. and Ford Global Technologies, Inc.

if nargin == 1
    TSModels= mdev.TwoStage;
else

    L= model(mdev);

    X = getdata(mdev,'FIT');
    XGd= double(X{end});
    prf= children(mdev);

    % get datum model
    DatumType= get(L,'DatumType');
    switch DatumType
        % initialise datum model for press predictions
        case {1,2}
            % datum is first child of local
            DatumModel= prf(1).model;
            Yrf= double(prf(1).getdata('Y'));
            DatumModel= iInitStore(DatumModel,XGd,Yrf);
        case 3
            % model is contained in datum link
            pdatum= datumlink(mdev);
            DatumModel= pdatum.model;
            Yrf= double(pdatum.getdata('Y'));
            DatumModel= iInitStore(DatumModel,XGd,Yrf);
        otherwise
            % no datum model
            DatumModel= 0;
    end

    % local RMSE
    sse=  mdev.MLE.SSE;
    df= mdev.MLE.df;
    DFT= sum(df(mdev.FitOK));
    if DFT
        s2= sum(sse(mdev.FitOK))/DFT;
    else
        s2= 0;
    end

    % response feature data
    YrfData= children(mdev,@getdata,'Y');
    nrf= length(YrfData);
    Yrf= zeros(size(YrfData{1},1), nrf);
    for i=1:nrf
        Yrf(:,i)= double(YrfData{i});
    end
    mrf= children(mdev,@model);
    % response feature start position 
    RF1= RFstart(L);
    for j= RF1+1:nrf;
        % initialise response feature models for press predictions
        mrf{j}= iInitStore(mrf{j},XGd,Yrf(:,j));
    end

    TSok= false( size(selrf,1),1 );
    TSModels= cell(size(TSok));
    for i=1:size(selrf,1);
        % select response features for localmod
        Lf= EvalDelG(SelFeat(L,selrf(i,:)));
        % select global models
        G= mrf(RF1+selrf(i,:));
        TS= xregtwostage(Lf,G,DatumModel);
        if canMLE(TS)
            % Check that a two-stage model can be built
            
            srf= selrf(i,:);
            % use local covariance component
            W= mdev.MLE.Sigma(srf,srf,:);
            % initialise PEV for univariate model
            [TSModels{i},TSok(i)]= pevinit(TS,XGd,Yrf(:,srf+RF1),s2*W,true);
        else
            TSModels{i} = TS;
            TSok(i) = canGTS(TS);
        end
    end
    
    TSModels= TSModels(TSok);
    selrf= selrf(TSok,:);
    mdev.TwoStage= TSModels;

    OldBest= bestmdev(mdev);
    if OldBest~=0
        OldRF=mdev.ResponseFeatures(OldBest,:);
        ind= ismember(selrf,OldRF,'rows');
        if any(ind)
            % Keep old best model - this will be selected in the model
            % selection tool
            mdev.modeldev= BestModel(mdev.modeldev,find(ind));
        else
            % No best model - allow model selection tool to suggest one
            mdev= BestModel(mdev,0);
        end
    end
    mdev.ResponseFeatures   = selrf;

    xregpointer(mdev);
end


function m= iInitStore(m,X,Yrf)

if islinear(m)
    % This is only needed for press predictions 
    X= double(X);
    Yrf= double(Yrf);
    % Need coded and transformed data
    [Xcode,Ytrans,DataOK]= checkdata(m,X,Yrf);
    % Need to set unused data to NaN so that press predictions will be
    % aligned with data
    Yrf(~DataOK)= NaN;
    X(DataOK,:)= Xcode;
    % Initstore without variance calculation
    m= InitStore(m,X,Yrf,~DataOK,false);
end