clear
reset(0)
echo on
%Solution to Problem 18

%In this problem we were asked to analyze
%some heat transfer data.  We copy and store
%the data under the file name heat.dat, and
%then load it in:
load heat.dat
re=heat(:,1);
pr=heat(:,2);
nu=heat(:,3);
pause
%We wish to fit the data to a power law
%model of the form:

%Nu = c * Re^a * Pr^b

%To use linear regression we must linearize
%the problem.  We do this by taking the log.
%Thus we have the regression problem:
n=length(nu);
m=3;
a=[log(re),log(pr),ones(n,1)];
b=log(nu);
%which has the solution:
x=a\b
%where x(1) is a, x(2) is b, and x(3) is log(c).
pause
%We may also calculate the matrix of covariance:
k=inv(a'*a)*a';
r=b-a*x;
varb=r'*r/(n-m)
varx=varb*k*k'
pause
%We thus have the values for x and the errors:
[x,diag(varx).^0.5]
pause
%Now we do it again using the bootstrap method.
%We shall take p samples of n data points from
%our replicated data set.  We do this using the
%random number generator provided by matlab:
p=500;
for i=1:p
 index=ceil(rand(n,1)*n);
 bboot=b(index);
 aboot=a(index,:);
 xboot(:,i)=aboot\bboot;
 echo off
end
echo on
pause
%We now have an array of solutions for our random
%samples of our replicated data set.  Let's look at
%the average of x:
avexboot=(sum(xboot')')/p
%which is close to the value which we had before:
ratio=avexboot./x
pause
%And then the matrix of covariance:
for i=1:m
 for j=1:m
  varxboot(i,j)=(xboot(i,:)*xboot(j,:)'/p- ...
  avexboot(i)*avexboot(j))*p/(p-1);
  echo off
 end
end
echo on
varxboot
pause
%Finally, we compare this to varx:
ratiovar=varxboot./varx
%which matches pretty closely.  Note that the
%value of p which we have chosen is much larger
%than n.  We need to have a large value for p
%because it is much harder to estimate the value
%of the variance than it is to estimate the mean
%(you need alot more points).  Remember that
%because the bootstrap uses a random sample, it
%(the means and variances) will change every
%time you run it!




__________________________




Save this data set under the file name heat.dat

   2.0000e+04   8.1600e+00   1.4318e+02
   2.0500e+04   8.3900e+00   1.6686e+02
   2.1000e+04   8.5600e+00   1.5335e+02
   2.1500e+04   7.1900e+00   1.5671e+02
   2.2000e+04   8.0500e+00   1.5221e+02
   2.2500e+04   7.6900e+00   1.5299e+02
   2.3000e+04   7.0500e+00   1.6264e+02
   2.3500e+04   7.4100e+00   1.7429e+02
   2.4000e+04   7.2400e+00   1.5987e+02
   2.4500e+04   7.0100e+00   1.5884e+02
   2.5000e+04   8.3700e+00   1.6875e+02
   2.5500e+04   7.1900e+00   1.7029e+02
   2.6000e+04   7.9500e+00   1.8549e+02
   2.6500e+04   8.3400e+00   1.8124e+02
   2.7000e+04   7.2400e+00   1.6957e+02
   2.7500e+04   7.4500e+00   1.8350e+02
   2.8000e+04   8.2800e+00   1.7738e+02
   2.8500e+04   8.5700e+00   2.0250e+02
   2.9000e+04   8.9400e+00   1.9098e+02
   2.9500e+04   7.8600e+00   1.8808e+02
   3.0000e+04   7.2600e+00   2.0125e+02
   3.0500e+04   7.2000e+00   1.9637e+02
   3.1000e+04   7.5800e+00   2.0122e+02
   3.1500e+04   7.0500e+00   1.9261e+02
   3.2000e+04   8.6600e+00   2.0215e+02
   3.2500e+04   8.8100e+00   2.1756e+02
   3.3000e+04   7.0200e+00   2.1018e+02
   3.3500e+04   8.1000e+00   2.2368e+02
   3.4000e+04   7.2300e+00   2.0038e+02
   3.4500e+04   8.2500e+00   2.1413e+02
   3.5000e+04   8.0700e+00   2.2018e+02
   3.5500e+04   8.9800e+00   2.3234e+02
   3.6000e+04   8.6200e+00   2.4733e+02
   3.6500e+04   7.9700e+00   2.3147e+02
   3.7000e+04   7.6100e+00   2.3418e+02
   3.7500e+04   7.6000e+00   2.3480e+02
   3.8000e+04   8.3700e+00   2.5935e+02
   3.8500e+04   8.9000e+00   2.4134e+02
   3.9000e+04   7.5100e+00   2.2648e+02
   3.9500e+04   7.7800e+00   2.4771e+02
   4.0000e+04   7.7700e+00   2.4293e+02