Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: Linear regression.
Date: Thu, 13 Mar 2008 18:50:19 +0000 (UTC)
Organization: Bedford Institute of Oceanography
Lines: 108
Message-ID: <frbt1b$9nn$1@fred.mathworks.com>
References: <fr9c9t$o1p$1@news.onet.pl>
Reply-To: <HIDDEN>
NNTP-Posting-Host: webapp-03-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1205434219 9975 172.30.248.38 (13 Mar 2008 18:50:19 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Thu, 13 Mar 2008 18:50:19 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1149396
Xref: news.mathworks.com comp.soft-sys.matlab:457110



Hello,

I had the same problem and did not find anything in matlab.
Below is the code, which I wrote up for this purpose. It 
contains the reference to the source, which I used. Perhaps 
I will consider uploading it to file exchange, but for now 
I will just copy-paste it. Hope you will cope with text 
wrapping and figure it out.

Yuri
------------------------
%*******************************
function [sesl, cisl, sein, ciin] = secisi(X,Y,pr);

% Returns Standard Errors and Confidence Intervals
% for the Slope and Intercept of a regression.

% Author: Yuri Geshelin, 2002.

% Source: N. R. Draper and H. Smith, 1966:
% Applied Regression Analysis, 407 p. Publ. by
% John Wiley and Sons, Inc. New York - London - Sydney.

% Input: 
% 1-D arrays X and Y, for which the regression
%        Y = B0 + B1 * X is computed
% pr - probability (e.g. if we want 95% 
%      confidence level, then pr = 95.

% Output: 
% sesl - standard error for slope B1
% clsl - pr%-confindence interval for slope B1
% sein - standard error for intercept B0 
% clin - pr%-confindence interval for intercept B0

if pr < 10 | pr > 99.9
   error('pr should be between 10 and 99.9!') 
end

n = length(X);
if length(Y) ~= n
   error('X, Y should be the same lengths!')
end

n2 = n - 2; % n - 2 degrees of freedom

if n2 < 121 % if n2 = 120, then we already would 
            % have to use the next-to-last line 
            % of the table on p. 305: this code is only
            % for the last line.
   error('This code only works with n > 120!')
end

% The following line is the result of 
% prob = 100 * ...
% (1 - [0.9:-0.2:0.3 0.2 0.1 0.05 0.02 0.01 0.001]) ...
% (see p. 305):

prob = [10 30 50 70 80 90 95 98 99 99.9]; % percentage

pp = [0.126 0.385 0.674 1.036 1.282 ...
      1.645 1.960 2.326 2.576 3.291]; % last line of the 
table,
                                      % p. 305.

% Compute total corrected sum of squares (SS) of the Y's
% (see p.14 and table on p. 16):

SStot = sum((Y - mean(Y)).^2);

% Compute SS due to regression (see p. 16):

SX2 = sum(X .^ 2);
SSrgr = ...
(sum(X .* Y) - sum(X) * sum(Y) / n) ^2 / ...
(SX2 - sum(X) ^2 / n);

% Compute SS due to residual (see p. 16):

SSres = SStot - SSrgr;

% Compute mean square about regression 
% (see tables on pp. 15, 16):

s2 = SSres / n2; 

% Compute estimate of variance of slope B1 (p. 19)

sesl = sqrt(s2 / sum((X - mean(X)) .^2 ));

% Found the interpolated value from the table
% (page 305):

val_int = interp1(prob,pp,pr);

% Compute confidence limit for slope:

cisl =  val_int * sesl;

% Compute estimate of variance of intercept B0 (p. 21)

sein = sesl * sqrt(SX2 / n);

% Compute confidence limit for intercept:

ciin = val_int * sein;
%***********************************