www.gusucode.com > mbctools 工具箱 matlab 源码程序 > mbctools/@mdev_local/Change_RespFeat.m
function mdev= Change_RespFeat(mdev,SweepNos) %CHANGE_RESPFEAT response feature data update % % mdev= Change_RespFeat(mdev); % change all sweeps/rfs % Copyright 2000-2015 The MathWorks, Inc. and Ford Global Technologies, Inc. if nargin<2 SweepNos = 1:size(mdev.AllModels,2); end L= model(mdev); % this is now a way to get clean data [X,Yresp]= getdata(mdev); [X,Y,DataOK,bd]= checkdata(L,X,Yresp); Yresp(bd)= NaN; % bad sweeps Yresp(:,:,~DataOK)= NaN; if any(DataOK) X(:,:)= invcode(L,double(X)); end if ~isempty(mdev.GLSWeights) Wc= mdev.GLSWeights; else Wc= cell(1,size(X,3)); end if isnull(mdev.RFData) % make response feature pointer if required XF =getdata(mdev,'fit'); XG = XF{2}; S= sweepset('Variable', 1 , '%g' , {}, 'Response Feature' , {} , 'none' , zeros(size(XG,1),0) ); D= [XG(:,~1) S]; pD=xregpointer(D); mdev.RFData = pD; end Nf= numfeats(L); Ns= length(SweepNos); NewRF=zeros(Ns,Nf); nl= nfactors(L); D= zeros(Ns,nl); if isempty(mdev.MLE) || size(mdev.AllModels,2) ~= size(mdev.MLE.Sigma,3) % size of data may have changed Ns= size(mdev.AllModels,2); Sigma = zeros(Nf,Nf,Ns); SSE_nat = zeros(Ns,1); sse = zeros(Ns,1); df = ones(Ns,1); ValRMSE = NaN(Ns,1); ValDataSize= zeros(Ns,1); else Sigma = mdev.MLE.Sigma; SSE_nat = mdev.MLE.SSE_nat; sse = mdev.MLE.SSE; df = mdev.MLE.df; if ~isfield(mdev.MLE,'ValRMSE') % initialise validation data storage mdev.MLE.ValRMSE = NaN(1,size(mdev.AllModels,2)); mdev.MLE.ValDataSize = NaN(1,size(mdev.AllModels,2)); end ValRMSE = mdev.MLE.ValRMSE; ValDataSize= mdev.MLE.ValDataSize; end if Nf~=size(Sigma,1) Sigma = zeros(Nf,Nf,Ns); end mdev.IsLinearised = false; m=1; [Xval,Yval]= getLocalValidationData(mdev,SweepNos); for testCounter= 1:length(SweepNos) testNumber = SweepNos(testCounter); if mdev.FitOK(testNumber) j = sum(DataOK(1:testNumber)); L= LocalModel(mdev,testNumber); Xc= code(L,X{j}); D(m,:)= datum(L); NewRF(m,:)= evalfeatures(L); % need to re-evaluate covariance matrix as well if ~isempty(Wc{testNumber}) && length(Y{j})~=size(Wc{testNumber},1) Wc{testNumber} = CalcWeights(L,Xc); mdev.GLSWeights{testNumber}= Wc{testNumber}; end [SummStats,sig]= localstats(L,Xc,Y{j},Wc{testNumber}); SSE_nat(testNumber) = SummStats.SSE_natural; sse(testNumber) = SummStats.sse; df(testNumber) = SummStats.df; if ~isempty(Xval) && ~isempty(Xval{testCounter}) % only calculate validate RMSE if available [ValRMSE(testNumber),ValDataSize(testNumber)] = LocalValidationRMSE(mdev,testNumber,L,Xval(testCounter),Yval(testCounter)); end if ~isempty(sig) Sigma(:,:,testNumber)= sig; end else SSE_nat(testNumber) = Inf; sse(testNumber) = Inf; df(testNumber) = 0; NewRF(m,:) = NaN; D(m,:) = NaN; end m=m+1; end mdev.MLE.Sigma = Sigma ; mdev.MLE.SSE_nat= SSE_nat; mdev.MLE.SSE = sse; mdev.MLE.df = df; mdev.MLE.ValRMSE = ValRMSE; mdev.MLE.ValDataSize = ValDataSize; if RFstart(L) NewRF= [D NewRF]; end % set up guids for response feature data and remove any global outliers which change sguids = getSweepGuids(Yresp,'goodonly'); RF = mdev.RFData.info; RF = setGuids(RF,sguids); % update rf data if ~isempty(NewRF) RF(SweepNos,:)= NewRF; end mdev.RFData.info= RF; mdev = updateOutlierIndices(mdev,Yresp,RF); % Update statistics NotFitted = ~mdev.FitOK; mdev.MLE.SSE(NotFitted)= 0; mdev.MLE.SSE_nat(NotFitted)= 0; mdev.MLE.df(NotFitted)= 0; mdev.MLE.ValRMSE(NotFitted) = NaN; mdev.MLE.ValDataSize(NotFitted) = 0; SSE= mdev.MLE.SSE; df = mdev.MLE.df; % Summary Statistics for local model DFT= sum(df(mdev.FitOK)); if DFT>0 Pooled_MSE= sum(SSE(mdev.FitOK))/DFT; else Pooled_MSE= NaN; end mdev.MLE.Pooled_MSE= Pooled_MSE; SSE= mdev.MLE.SSE_nat; ind= isfinite(SSE(:)) & mdev.FitOK(:); DFT= sum(df(ind)); if DFT>0 Pooled_MSE= sum(SSE(ind))/DFT; else Pooled_MSE= NaN; end % mdev_local summary stats ch = colhead(mdev); stats = [sqrt(Pooled_MSE) NaN(1,length(ch)-1)]; if isPointByPoint(mdev) % point-by-point validation statistics sse = mdev.MLE.ValRMSE.^2.*mdev.MLE.ValDataSize; HasVal = isfinite(sse); if any(HasVal) ValRMSE = sqrt(sum(sse(HasVal))/sum(mdev.MLE.ValDataSize(HasVal))); stats(2) = ValRMSE; end end mdev= statistics(mdev,stats); pointer(mdev);