I suggest to use non-scalar structures:
% Some random inputs
y = rand(501,1);
sincos = rand(5,1);
n = 10;
% Preallocate structure
data(1:5) = struct('X',NaN(501,3),'bhat',NaN(3,1),'yhat',NaN);
% LOOP
for ii = 1:5
tmp = 2*pi*(sincos(ii))*n;
data(ii).X = rand(501,3);
data(ii).X(:,2) = cos(tmp)';
data(ii).X(:,3) = sin(tmp)';
data(ii).bhat = data(ii).X\y;
data(ii).yhat = data(ii).bhat(1)+ data(ii).bhat(2)*cos(tmp) + data(ii).bhat(3)*sin(tmp);
end
Now, each element of the structure array has the same fields but correspond to a different regression. You don't need to hardcode or play around with different names, just remember that the 5th element of the structure has the the 5th regression.
Also, another advantage of this way of organizing data is that you can "collect" the data in the same fields in a very easy way:
[data.bhat]
will concatenate horizontally the estimated coefficients. Each column represents a different regression. With the preallocation you can mix different models (with different number of regressors) and still be able to concatenate results (you can also place them in the right position by hardcoding a vector with positions).
