www.pudn.com > NumericalComputing.rar > censusgui.m, change:2003-10-25,size:6530b


function censusgui(callbackarg) 
%CENSUSGUI Try to predict the US population in the year 2010. 
% This example is older than MATLAB.  It started as an exercise in 
% "Computer Methods for Mathematical Computations", by Forsythe, 
% Malcolm and Moler, published by Prentice-Hall in 1977. 
% The data set has been updated every ten years since then. 
% Today, MATLAB makes it easier to vary the parameters and see the 
% results, but the underlying mathematical principles are unchanged: 
% 
%    Using polynomials of even modest degree to predict 
%    the future by extrapolating data is a risky business. 
% 
% The data is from the decennial census of the United States for the 
% years 1900 to 2000.  The task is to extrapolate beyond 2000. 
% In addition to polynomials of various degrees, you can choose 
% interpolation by a cubic spline, interpolation by a shape-preserving 
% Hermite cubic, and a least squares fit by an exponential. 
% Error estimates attempt to account for errors in the data, 
% but not in the extrapolation model. 
 
% Census data for 1900 to 2000. 
% The population on April 1, 2000 was 281,421,906, according to: 
% http://www.census.gov/main/www/cen2000.html 
 
p = [ 75.995  91.972 105.711 123.203 131.669 150.697 ... 
     179.323 203.212 226.505 249.633 281.422]'; 
 
t = (1900:10:2000)';   % Census years 
x = (1890:1:2019)';    % Evaluation years 
w = 2010;              % Extrapolation target 
guess = 320;           % Eyeball extrapolation 
z = guess;             % Extrapolated value 
dmax = length(t)-1;    % Maximum polynomial degree 
 
if nargin == 0 
 
   % Initialize plot and uicontrols 
 
   shg 
   clf reset 
   set(gcf,'doublebuffer','on','name','Census gui', ... 
       'menu','none','numbertitle','off') 
   h.plot = plot(t,p,'bo', x,0*x,'k-', w,0,'.', [x;NaN;x],[x;NaN;x],'m:'); 
   darkgreen = [0 2/3 0]; 
   darkmagenta = [2/3 0 2/3]; 
   marksize = get(0,'defaultlinemarkersize'); 
   set(h.plot(3),'color',darkgreen,'markersize',4*marksize-6) 
   set(h.plot(4),'color',darkmagenta) 
   axis([min(x) max(x) 0 400]) 
   title('Predict U.S. Population in 2010') 
   ylabel('Millions') 
 
   h.text = text(w-16,z+10,'predict','color',darkgreen,'fontweight','bold'); 
   h.model = uicontrol('units','norm','pos',[.20 .80 .20 .05], ... 
           'style','popup','background','white','string', ... 
           {'census data','polynomial','pchip','spline','exponential'}, ... 
           'callback','censusgui([])'); 
   h.deg = uicontrol('units','norm','pos',[.26 .75 .13 .04], ... 
           'tag','degree','style','text','background','white', ... 
           'userdata',3,'string','degree = 3'); 
   h.ls = uicontrol('units','norm','pos',[.20 .75 .05 .04], ... 
           'style','push','string','<','fontweight','bold', ... 
           'callback','censusgui(''<'')'); 
   h.gt = uicontrol('units','norm','pos',[.40 .75 .05 .04], ... 
           'style','push','string','>','fontweight','bold', ... 
           'callback','censusgui(''>'')'); 
   h.err = uicontrol('units','norm','pos',[.20 .65 .20 .05], ... 
           'style','check','background','white', ... 
           'string','error estimates','callback','censusgui([])'); 
   set(gcf,'userdata',h); 
   uicontrol('style','push','units','normal','pos',[.85 .02 .10 .06], ... 
      'string','close','callback','close(gcf)') 
   callbackarg = []; 
 
else 
 
   h = get(gcf,'userdata'); 
  
end 
 
% Polynomial degree 
 
d = get(h.deg,'userdata'); 
if isequal(callbackarg,'<'), d = d - 1; end 
if isequal(callbackarg,'>'), d = d + 1; end 
set(h.deg,'userdata',d) 
 
% Update plot with new model 
 
models = get(h.model,'string'); 
model = models{get(h.model,'value')}; 
switch model 
   case 'census data' 
      y = NaN*x; 
      z = 320; 
   case 'polynomial' 
      s = (t-1950)/50;   c = polyfit(s,p,d); 
      s = (x-1950)/50;   y = polyval(c,s); 
      s = (w-1950)/50;   z = polyval(c,s); 
   case 'pchip' 
      y = pchip(t,p,x); 
      z = pchip(t,p,w); 
   case 'spline' 
      y = spline(t,p,x); 
      z = spline(t,p,w); 
   case 'exponential' 
      c = polyfit(log(t),log(p),1); 
      y = exp(polyval(c,log(x))); 
      z = exp(polyval(c,log(w))); 
end 
set(h.plot(2),'ydata',y); 
set(h.plot(3),'ydata',z); 
set(h.text,'pos',[w-18,min(max(z+10,20),380)],'string',sprintf('%8.3f',z)) 
 
% Update controls 
 
switch model 
   case 'census data' 
      set(h.err,'vis','off','value',0); 
      set([h.deg; h.gt; h.ls],'vis','off'); 
      set(h.text,'pos',[w-16,z+10],'string','predict') 
   case 'polynomial' 
      set(h.err,'vis','on','pos',[.20 .68 .20 .05]); 
      set(h.deg,'vis','on','string',['degree = ' num2str(d)]); 
      set([h.gt; h.ls],'vis','on','enable','on'); 
      if d == 0, set(h.ls,'enable','off'), end 
      if d == dmax, set(h.gt,'enable','off'), end 
   otherwise 
      set(h.err,'vis','on','pos',[.20 .75 .20 .05]); 
      set([h.deg; h.gt; h.ls],'vis','off'); 
end 
 
% Display error estimates if requested 
 
if get(h.err,'value') == 1 
   errest = errorestimates(model,t,p,x,y,d); 
   set(h.plot(4),'vis','on','ydata',errest); 
else 
   set(h.plot(4),'vis','off'); 
end 
 
 
% ------------------------------------------------ 
 
function errest = errorestimates(model,t,p,x,y,d) 
% Provide error estimates for censusgui 
 
switch model 
   case 'polynomial' 
      if d > 0 
         V(:,d+1) = ones(size(t)); 
         s = (t-1950)/50; 
         for j = d:-1:1 
            V(:,j) = s.*V(:,j+1); 
         end 
         [Q,R] = qr(V); 
         R = R(1:d+1,:); 
         RI = inv(R); 
         E = zeros(length(x),d+1); 
         s = (x-1950)/50; 
         for j = 1:d+1 
            E(:,j) = polyval(RI(:,j),s); 
         end 
         sig = 10;   % Rough estimate 
         e = sig*sqrt(1+diag(E*E')); 
         errest = [y-e; NaN; y+e]; 
      else 
         errest = [y-NaN; NaN; y+NaN]; 
      end 
   case {'pchip','spline'} 
      n = length(t); 
      I = eye(n,n); 
      E = zeros(length(x),n); 
      for j = 1:n 
         if isequal(model,'pchip') 
            E(:,j) = pchip(t,I(:,j),x); 
         else 
            E(:,j) = spline(t,I(:,j),x); 
         end 
      end 
      sig = 10;  % Rough estimate 
      e = sig*sqrt(1+diag(E*E')); 
      errest = [y-e; NaN; y+e]; 
   case 'exponential' 
      V = [ones(size(t)) log(t)]; 
      [Q,R] = qr(V); 
      c = R\(Q'*log(p)); 
      r = log(p) - V*c; 
      E = [ones(size(x)) log(x)]/R(1:2,1:2); 
      sig = norm(r); 
      e = sig*sqrt(1+diag(E*E')); 
      errest = [y.*exp(-e); NaN; y.*exp(e)]; 
end