While robustfit gives as an output robust weights, it does not have an option of using a custom fit function as an input.
Since I needed to fit specific models to my data, I used instead fit, which is more flexible int hat sense, with robust option activated. However, fit does not print out robust fit related statistics, nor weights. I would like to have access to weights, in order to compare several models via weighted least square error. I tried to calculate weights myself, but the outcome does not seem right.
ft = fittype(......);
opts = fitoptions(......);
opts.Robust = 'Bisquare';
[fitresult, gof, output] = fit(x,y,ft,opts);
res = output.residuals;
[Q,~] = qr(output.Jacobian,0);
h = min(.9999, sum(Q.*Q,2)); %leverage
res_adj = res./sqrt(1-h); %adjusted residuals
K = 4.685; %tuning constant
res_adj_sorted = sort(abs(res_adj));
s = median(res_adj_sorted(max(1,p):end)) / 0.6745;
u=res_adj./(K*s); % standardized residuals
w = zeros(size(res));
for i = 1:length(res)
if abs(u(i,:)) < 1
w(i,:) = (1-u(i,:)^2)^2;
else
w(i,:) = 0;
end
end
