0

I've been running a variation of the doseResponse function downloaded from here to generate dose-response sigmoid curves. However, I've had trouble with one of my datasets generating a linear curve instead. By running the following code, I get the following error and produce the following graph. I also uploaded the data called dose2.csv and resp2.csv to google drive here. Does anyone know how I can fix this? Thanks.

Code to generate graph

% Plotting Dose-Response Curve
response = resp2;
dose = dose2; 

% Deal with 0 dosage by using it to normalise the results.
normalised=0;
if (sum(dose(:)==0)>0)
    %compute mean control response
    controlResponse=mean(response(dose==0));
    %remove controls from dose/response curve
    response=response(dose~=0)/controlResponse;
    dose=dose(dose~=0);
    normalised=1;
end

%hill equation sigmoid
sigmoid=@(beta,x)beta(1)+(beta(2)-beta(1))./(1+(x/beta(3)).^beta(4));

%calculate some rough guesses for initial parameters
minResponse=min(response);
maxResponse=max(response);
midResponse=mean([minResponse maxResponse]);
minDose=min(dose);
maxDose=max(dose);

%fit the curve and compute the values
%[coeffs,r,J]=nlinfit(dose,response,sigmoid,[minResponse maxResponse midResponse 1]); % nlinfit doesn't work as well
beta_new = lsqcurvefit(sigmoid,[minResponse maxResponse midResponse 1],dose,response);
[coeffs,r,J]=nlinfit(dose,response,sigmoid, beta_new);

ec50=coeffs(3);
hillCoeff=coeffs(4);

%plot the fitted sigmoid
xpoints=logspace(log10(minDose),log10(maxDose),1000);
semilogx(xpoints,sigmoid(coeffs,xpoints),'Color',[1 0 0],'LineWidth',2)
hold on

%notate the EC50
text(ec50,mean([coeffs(1) coeffs(2)]),[' \leftarrow ' sprintf('EC_{50}=%0.2g',ec50)],'FontSize',20,'Color',[1 0 0]);

%plot mean response for each dose with standard error
doses=unique(dose);
meanResponse=zeros(1,length(doses));
stdErrResponse=zeros(1,length(doses));
for i=1:length(doses)
    responses=response(dose==doses(i));
    meanResponse(i)=mean(responses);
    stdErrResponse(i)=std(responses)/sqrt(length(responses));
    %stdErrResponse(i)=std(responses);
end

errorbar(doses,meanResponse,stdErrResponse,'o','Color',[1 0 0],'LineWidth',2,'MarkerSize',12)

Warning Message

Solver stopped prematurely.

lsqcurvefit stopped because it exceeded the function evaluation limit,
options.MaxFunctionEvaluations = 4.000000e+02.

Warning: Iteration limit exceeded.  Returning results from final iteration. 

Graph (looking to generate a sigmoid curve not linear) enter image description here

won5830
  • 512
  • 5
  • 13

1 Answers1

0

You also need to optimize your initial value [minResponse maxResponse midResponse 1] for lsqcurvefit. Don't just simply start with minimum or maximum values of given values. Instead, you may first start with your equations to estimate your coefficients.

Given the sigmoid model of sigmoid=@(beta,x)beta(1)+(beta(2)-beta(1))./(1+(x/beta(3)).^beta(4)). As x gets arbitrarily close to inf, equation will return beta(2). And as x gets arbitrarily close to 0, equation will return beta(1). Therefore, initial estimation of minResponse, maxResponse, and midResponse seems reasonable enough. Actually your problem lies in your initial estimation of 1. beta(4) can be roughly estimated with the inclination of your log graph. To my rough sketch it was around 1/4 and therefore you may conclude that your initial estimation of 1 was too large for convergence.

beta_new = lsqcurvefit(sigmoid,[minResponse maxResponse midResponse 1/4],dose,response);

Revised

won5830
  • 512
  • 5
  • 13