Fitting Weibull function using nlinfit - errors for some data not others

Hi there,

I am trying to fit a curve in MATLAB R2020B to some human perception data (psychometric function). The xaxis has different difficulties of the perception task (dataStims) and the yaxis has the proportion of correct responses (dataProp), where chance = .5.

Data examples from two participants are below

Dataset 1

dataStims = [10, 14, 25, 40, 55, 66, 70];

dataProp = [0.94444, 1, 0.77778, 0.88889, 0.55556, 0.44444, 0.66667];

Dataset 2

dataStims = [10, 14, 25, 40, 55, 66, 70]; % these are the same for all participants

dataProp = [1, 1, 0.777777777777778, 0.722222222222222, 0.722222222222222, 0.500000000000000, 0.555555555555556];

I am currently trying to fit a Weibull function to this data. Function below

% Weibull with two free parameters (gamma and lambda are set - add their values below)

gamma = 0.5; lambda = 0.01;

f = @(F,x) gamma + (1 - gamma - lambda).*(1 - exp(-1*(x./F(1)).^F(2)));

Using nlinfit to do the fitting. Code below

startvalues=[1 1];

options = statset('nlinfit');

options.RobustWgtFun = [ ];% 'bisquare';

[params, ~] = nlinfit(dataStims,dataProp,f,startvalues,options,'Weights',dataCount);

For about half of my participants (e.g., dataset 1) this works well, but for the other half it fails (e.g., dataset 2) and gives this error:

Error using nlinfit>checkFunVals (line 649)

The function you provided as the MODELFUN input has returned Inf or NaN values.

Error in nlinfit>LMfit (line 596)

if funValCheck && ~isfinite(sse), checkFunVals(r); end

Error in nlinfit (line 284)

[beta,J,~,cause,fullr] = LMfit(X,yw, modelw,beta,options,verbose,maxiter);

Error in s3_fit_Weibulls (line 28) [params, ~] = nlinfit(dataStims,dataProp,f,startvalues,options,'Weights',dataCount);

Error in s0_run_fitting_240222 (line 99) s3_fit_Weibulls

I don't know enough about fitting functions or nlinfit to properly diagnose what is happening and was wondering if anyone could help?

I do know, obviously, that the data is noisy - this is just data from a quick pilot, the real data will have way more trials. I presume the Weibull function is just not very tolerant to noisy data, and this is why it fails? Is there a way to make it more robust to noise and run on more datasets, or something similar? I saw in the nlinfit documentation that you can change some of the control parameters used by nlinfit

% 'MaxIter' - Maximum number of iterations allowed. Defaults to 100.

% 'TolFun' - Termination tolerance on the residual sum of squares.

% Defaults to 1e-8.

% 'TolX' - Termination tolerance on the estimated coefficients

% BETA. Defaults to 1e-8.

Does anyone know whether I should try changing any of these and what it would mean? If I should try it, I can't work out how to put the options into the call for nlinfit

I tried this code, but it does not work:

[params, ~] = nlinfit(dataStims,dataProp,f,startvalues,options,'MaxIter', 200, 'Weights',dataCount);

Alternately, if anyone has a suggestion for a function/ curve that might better fit data that looks like mine (assumptions: y values, chance = 0.5, max performance = 1.0, approximately S shaped decreasing function), then that would be very much appreciated.

I first tried fitting a logistic regression, but it does not seem to fit as my data spans from ~0.5 on the yaxis to 1.0, and it seems to want 0-1.0 range to make a S shaped curve that fits the data well.

Any help would be very very greatly appreciated!

Thank you in advance,

Harriet

Answer 1

[Torsten]

Start with

startvalues=[35 1];

instead of

startvalues=[1 1];

Further, I suggest you choose

f = @(F,x) (1 - exp(-(x./F(1)).^F(2)));

instead of

f = @(F,x) gamma + (1 - gamma - lambda).*(1 - exp(-1*(x./F(1)).^F(2)));

Here is a different solver you could test:

dataStims = [10, 14, 25, 40, 55, 66, 70]; % these are the same for all participants
dataProp = [1, 1, 0.777777777777778, 0.722222222222222, 0.722222222222222, 0.500000000000000, 0.555555555555556];
%dataProp = [0.94444, 1, 0.77778, 0.88889, 0.55556, 0.44444, 0.66667];
%f = @(F) (1 - exp(-1*(dataStims./F(1)).^F(2))) - dataProp;
gamma = 0.5; lambda = 0.01;
f = @(F) gamma + (1 - gamma - lambda).*(1 - exp(-1*(dataStims./F(1)).^F(2))) - dataProp;
p0 = [35 1];
p = lsqnonlin(f,p0);
plot(dataStims,f(p)+dataProp)
hold on
plot(dataStims,dataProp,'+')

[HDJ]

Hi Torsten,

Thank you so much for your help - simply putting more appropriate starting values allowed all the fits to converge and produce a fit. Exactly as I was hoping for.

I just had a couple of questions to follow up some other parts of your answer:

Could you comment on why you suggest to use this function (f = @(F,x) (1 - exp(-(x./F(1)).^F(2)));) instead of the one I originally posted? The one I posted is for a 4-parameter Weibull, which allows you to set gamma and lambda and it estimates alpha and beta. Is the one you posted just the two parameter Weibull?

Also, what is the difference with the third 'different solver' method you posted at the bottom (f = @(F) gamma + (1 - gamma - lambda).*(1 - exp(-1*(dataStims./F(1)).^F(2))) - dataProp;). Could you comment on why this one is preferable?

Thank you for your help,

Harriet

[Torsten]

Now that I read your posting again, I get the feeling that fitting with the Weilbull function is not appropriate in your case. The usual 2-parameter Weilbull gives 0 at data_stims = 0 while you need a value of 1. I suggest

a*exp(-bx^c), e.g.

"lsqcurvefit" or "lsqnonlin" are other parameter fitting tools from the optimization toolbox I often use. I don't know exactly the advantages and disadvanteges compared to nlinfit.

Like Star Strider, I do not quite understand your approach to evaluate the data. What do you aim at with fitting the results of your experiment for each test person separately ?

[HDJ]

Hi Torsen,

Thanks for your response

What is this function? a*exp(-bx^c) - Could you link to some documentation about it?

In perceptual psychophysics, it is typical to fit a function to each individual participant's data. This allows you to get a curve that describes their performance on all stimuli. Then you extract a threshold to describe their performance at a certain point on the curve, e.g., what is the stimulus that corresponds with performance at (arbirarily) 50%/ 75%/ 82%