R = 0.3;
IT = 600:600:3000;
P = 0.97;
DL = R*IT.^P + 0.05*randn(size(IT));
% generate matrices:
pMAT = 0.01*randn(3)+P;
rMAT = 0.1*randn(3)+R;
dlMAT = bsxfun(@times,rMAT,bsxfun(@power,permute(IT,[3,1,2]),pMAT));
sol = cell(size(rMAT));
% regression:
for ind1 = 1:size(dlMAT,1)
for ind2 = 1:size(dlMAT,2)
sol{ind1,ind2} = polyfit(log(IT(:)),log(squeeze(dlMAT(ind1,ind2,:))),1);
end
end
fittedP = cellfun(@(x)x(1),sol);
fittedR = cellfun(@(x)exp(x(2)),sol);