function [ci,varb,corrb,varinf] = nlparci(x,f,J)
%NLPARCI Confidence intervals and Var/Corr on parameters of nonlinear models.
%   [CI,VARB,CORRB,VARINF] = NLPARCI(X,F,J) returns the 95% confidence interval CI
%   and the variance and correlation matrix for the parameters b
%   and VARINF, the variance inflation factors(>10 implies multicollinearity)
%   on the nonlinear least squares parameter estimate X, given the 
%   residual sum of squares, F, and the Jacobian matrix ,J, at the solution.
%   The confidence interval calculation is valid for systems where 
%   the number of rows of J exceeds the length of X. 
%
%   NLPARCI uses the outputs of NLINFIT for its inputs.
%   Example:
%      [x,f,J]=nlinfit(input,output,model,xinit);
%      [ci,varb,corrb] = nlparci(x,f,J);
%
%   See also NLINFIT.
%

%   Bradley Jones 1-28-94
%   Copyright (c) 1993-96 by the MathWorks, Inc.
%   $Revision: 2.3 $  $Date: 1996/02/13 23:26:30 $
%   Modified by aj June '96 to provide Var and Corr and Varinf info also

%Arthur Jutan 
%Dept of Chemical and Biochemical Enginering
%Univ of Western Ontario
%London Ontario Canada, N6A 5B9

%initialization
if nargin < 3
   error('Requires three inputs.');
end;

f = f(:);
[m,n] = size(J);
if m <= n
   error('The number of observations must exceed the number of parameters.');
end;

if length(x) ~= n
   error('The length of x must equal the number of colums in J.')
end

% approximation when a column is zero vector
temp = find(max(abs(j)) == 0);
if ~isempty(temp)
   J(temp,:) = J(temp,:) + sqrt(eps);
end;

%calculate covariance
[Q R] = qr(J,0);
Rinv = R\eye(size(R));
diag_info = sum((Rinv.*Rinv)')';

v = m-n;
rmse = sqrt(sum(f.*f)/v);

% calculate confidence interval
delta = sqrt(diag_info) .* rmse*tinv(0.975,v);
ci = [(x(:) - delta) (x(:) + delta)];

% Var and Corr Calculation for parameters --aj May '96
xtxinv=Rinv*Rinv';% inv(X'X) see regstat
varb=xtxinv*rmse*rmse;
%calc correlation matrix corrb
	varinf = diag(xtxinv); % variance inflation factors if >10 ==> collinearity
	corrb = xtxinv./sqrt(varinf*varinf');

%--end of nlparci.m---

