How do I fit negative values to a log-logistic probability distribution function?

To calculate the Standardised Precipitation Evapotranspiration Index (SPEI), a simple water budget(deficit, D) is calculated as the difference between precip (P) and potential evapotranspiration (ET), D=P-Et; and the deficits fitted to a log-logistic pdf to take care of the negative values, common in D the series, when ET is more than precipitation (Vicente, et al., 2010). http://journals.ametsoc.org/doi/pdf/10.1175/2009JCLI2909.1

I am writing a code in Matlab to calculate the SPEI. I have computed the deficit but I cannot fit negative deficits to a log-logistic pdf.

Am getting this error:

Error using prob.LoglogisticDistribution>loglfit (line 319) The data in X must be positive

Error in prob.LoglogisticDistribution.fit (line 165) p = loglfit(x,0.05,cens,freq,opt);

Any suggestion/ideas on how to go about it? Please help a brother.

Answer 1

[dpb]

The Matlab log-logistic distribution object isn't the extended three-parameter log-logistic distribution used in the subject paper. The author references a procedure for estimating the coefficients:

"Parameters of the log-logistic distribution can be obtained following different procedures. Among them, the L-moment procedure is the most robust and easy approach (Ahmad et al. 1988)."

I've never used the particular distribution in anger so can't say specifically what issues you might run into just trying the direct nonlinear least squares estimators built into Matlab, but the distribution fitted by the function you tried isn't intended for the problem you're trying to solve.

[Muiruri]

dpb, thank you for your response.

Stagge et al., (2015) suggest the GEV as the alternative (https://www.eea.europa.eu/data-and-maps/indicators/river-flow-drought-2/stagge-et-al-2015-future) though there's a rejoinder by Vicente-Serrano, et al., (2015). (<http://onlinelibrary.wiley.com/doi/10.1002/joc.4474/abstract>)

The GEV somewhat works well for my program and am getting comparable results.

SPEI was formulated in R originally. Thus, Matlab resources on SPEI are limited; although I still am convinced there's a simpler way to fit the L-moments and MLEs in Matlab. Maybe you have an idea or a script I can adopt?

Thanks in advance.