Linear Regression and Curve Fitting

4 ビュー (過去 30 日間)
Sam
Sam 2013 年 12 月 4 日
コメント済み: Wayne King 2013 年 12 月 4 日
I have a model and some data I'd like to fit to it: X_t = B1*cos(2*pi*omega*t) + B2*sin(2*pi*omega*t) + eta_t
What function would I use to conduct linear regression here, to find B1 and B2?

採用された回答

Wayne King
Wayne King 2013 年 12 月 4 日
編集済み: Wayne King 2013 年 12 月 4 日
Fs = 1000;
t = 0:1/Fs:1-1/Fs;
y = 1.5*cos(2*pi*100*t)+0.5*sin(2*pi*100*t)+randn(size(t));
y = y(:);
X = ones(length(y),3);
X(:,2) = cos(2*pi*100*t)';
X(:,3) = sin(2*pi*100*t)';
beta = X\y;
beta(1) is the estimate of the constant term, beta(2) the estimate of B1 and beta(3) the estimate of B2.
If you set the random number generator to its default for reproducible results:
rng default
Fs = 1000;
t = 0:1/Fs:1-1/Fs;
y = 1.5*cos(2*pi*100*t)+0.5*sin(2*pi*100*t)+randn(size(t));
y = y(:);
X = ones(length(y),3);
X(:,2) = cos(2*pi*100*t)';
X(:,3) = sin(2*pi*100*t)';
beta = X\y;
The results are:
beta =
-0.0326
1.5284
0.4643
pretty good.
  2 件のコメント
Sam
Sam 2013 年 12 月 4 日
Thank-you Wayne. In my case, where I have the time series data in a .dat file, am I right in thinking that I would load that in as 'y' in place of the y = "1.5*cos(2*pi*100*t)+0.5*sin(2*pi*100*t)+randn(size(t));"?
Wayne King
Wayne King 2013 年 12 月 4 日
Yes, that's correct. Make sure you flip it to a column vector if it isn't already. Obviously you have to change the frequencies in the design matrix to suit your problem.
If you are trying to estimate other frequencies you need to add two columns to the design matrix per frequency --- one for cosine and one for sine

サインインしてコメントする。

その他の回答 (1 件)

Jos (10584)
Jos (10584) 2013 年 12 月 4 日
You might be interested in the function REGRESS
X = [cos(2*pi*omega*t(:)) sin(2*pi*omega*t(:))]
B = regress(Y,X)
If you want to specify an offset B(3), add a column of ones to X
X(:,3) = 1

カテゴリ

Help Center および File ExchangeLinear and Nonlinear Regression についてさらに検索

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by