How to vectorize this loop ?

4 vues (au cours des 30 derniers jours)
Omer
Omer le 26 Juil 2012
I have a time series of X and Y variables and I want to run a linear regression on it based on a specific time window length. This ends up giving me a time series of regression coefficients. My current implmentation using a loop to picking out data and running a regression on it. How can I vectorize this piece of code?
window_length = 10;
regression_coeff = [];
for i =1:length(X) - window_length
regression_coeff = [regression_coeff; regress(Y(i:i+window_length ),X(i:i+window_length )];
end

Réponse acceptée

Teja Muppirala
Teja Muppirala le 27 Juil 2012
You are calling regress with two vectors. That is the same as a projection using the dot product. So then your code could be written like this:
window_length = 10;
W =ones(window_length+1,1);
regression_coeff = conv(X.*Y,W,'valid')./conv(X.^2,W,'valid');

Plus de réponses (1)

George Papazafeiropoulos
George Papazafeiropoulos le 27 Juil 2012
Modifié(e) : George Papazafeiropoulos le 29 Juil 2012
window_length = 10;
X(X==0)=Inf;
X_2=repmat(X,1,length(X)-window_length);
X_3=triu(X_2,-window_length)-triu(X_2,1);
X_4=reshape(X_3,length(X)*(length(X)-window_length),1);
X_4(X_4==0)=[];
X_4(X_4==Inf)=0;
X_5=X_4(~isnan(X_4));
Xfinal=reshape(X_5,1+window_length,length(X_5)/(1+window_length));
Y(Y==0)=Inf;
Y_2=repmat(Y,1,length(Y)-window_length);
Y_3=triu(Y_2,-window_length)-triu(Y_2,1);
Y_4=reshape(Y_3,length(Y)*(length(Y)-window_length),1);
Y_4(Y_4==0)=[];
Y_4(Y_4==Inf)=0;
Y_5=Y_4(~isnan(Y_4));
Yfinal=reshape(Y_5,1+window_length,length(Y_5)/(1+window_length));
regression_coeff = diag((Xfinal./repmat(sum(Xfinal.^2,1).^(1/2),1+window_length,1))'*Yfinal)./(sum(Xfinal.^2,1).^(1/2))';
Best regards,
George Papazafeiropoulos
  2 commentaires
Jan
Jan le 27 Juil 2012
Please, George, add the links to your profile, but not in your answers.
George Papazafeiropoulos
George Papazafeiropoulos le 29 Juil 2012
OK Jan, I removed them. Sorry for the inconvenience.

Connectez-vous pour commenter.

Produits

Community Treasure Hunt

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

Start Hunting!

Translated by