CODEWITHSUNDEEP
matlabregressionleast-squares
Matlabber
Matlabber

Reputation: 153

MATLAB implementation for Flexible Least Squares (FLS) regression

Is there a MATLAB implementation for Flexible Least Squares (FLS) regression1?

I am looking for a decent (well, the least painful) way to handle regression with time-varying coefficients. All ideas are welcome, but the best would be a pointer to a FLS regression implementation, because I couldn't find one after extensive googling.

1. R. Kalaba and L. Tesfatsion (1989), "Time-Varying Linear Regression via Flexible Least Squares", Computers and Mathematics with Applications, Vol. 17, pp. 1215-1245

Upvotes: 2

Views: 1856

Answers (1)

Yamaneko
Yamaneko

Reputation: 3563

I've found a R code on GitHub of the FLS regression. I have converted it to MATLAB, but it may have some errors due to differences between the functions of R and MATLAB. Now, we just need to test it in order to do the final adjustments.

function X = fls(A, b, mu, ncap, smoothed)

if isempty(mu)
    mu = 1;
end

if isempty(ncap)
    ncap = length(b);
end

if isempty(smoothed)
    smoothed = true;
end

[m n] = size(A);

M = zeros(n, n, ncap);
E = zeros(n, ncap);
X = zeros(n, ncap);
R = eye(n) * mu;

for j = 1:ncap
    Z = linsolve( qr( R + tcrossprod(A(j,:)) ), eye(n) );
    M(:,:,j) = mu * Z;               % (5.7b)
    v = b(j) * A(j,:);

    if (j == 1) 
        p = zeros(n,1);
    else
        p = mu * E(:,j-1);
    end

    w = p + v;

    E(:,j) = Z * w;                  % (5.7c)
    R = -mu * mu * Z;

    diag = R(1:n+1:end);
    R(diag) = R(diag) + 2*mu;
end

% Calculate eqn (5.15) FLS estimate at ncap
Q = -mu * M(:,:,ncap-1);

diag = Q(1: (size(Q,1)+1) : end);
Q(diag) = Q(diag) + mu;
Ancap = A(ncap,:);

C = Q + Ancap' * Ancap;
d = mu * E(:,ncap-1) + b(ncap)*Ancap';
X(:,ncap) = C * d;
X(:,ncap) = linsolve(qr(C),d);

if smoothed
    % Use eqn (5.16) to obtain smoothed FLS estimates for
    % X[,1], X[,2], ..., X[,ncap-1]

    for j = 1:ncap-1
        l = ncap - j;
        X(:,l) = E(:,l) + M(:,:,l) * X(:,l+1);
    end
else
    X = X(:,ncap);
end

end

function resp = tcrossprod(A)
resp = A * A';
end

Upvotes: 1

Related Questions