How do you derive (or where can I find a derivation of) the equations to match the five parameters of a general Pearson distribution to a mean, variance, skewness, kurtosis, and condition that the function integrate to one?

I have been using the method described in the Mathematica documentation under "Examples" and then "Applications". It provides the equations but no derivation. The four equations for the moments are set up as:

eq[r_] := r*Subscript[b, 0]*Moment[r - 1] + (r + 1) Subscript[b, 1]*Moment[r] 
+ (r + 2) Subscript[b, 2]*Moment[r + 1] - Subscript[a, 1]*Moment[r + 1] 
- Subscript[a, 0]*Moment[r]


meq = Table[MomentConvert[eq[r], CentralMoment], {r, 0, 3}] /. 
    {Moment[1] -> μ, CentralMoment[2] -> σ^2, CentralMoment[3] -> Sqrt[Subscript[β, 1]] σ^3, 
        CentralMoment[4] -> Subscript[β, 2] σ^4}

These set up four linear equations in five unknowns.

I am interested in the equation added to fix the normalizing coefficient (so the distribution integrates to one):

meq = Join[meq, {Subscript[a, 0] + (12 μ Subscript[β, 1] + 
   2 μ (9 - 5 Subscript[β, 2]) + σ Sqrt[Subscript[β, 1]] (3 + Subscript[β, 2]))}];

The equations for the moments themselves are pretty easy to derive, but I have not figured out the fifth equation for normalization.


1 Answer 1


We solve this in Chapter 5 (Section 5.2 E: Higher Order Pearson-Style Families) of Mathematical Statistics with Mathematica. To the best of my knowledge, the 5 parameter Pearson had not been solved prior to the first edition of our book (if someone is aware of an earlier solution, I'd be pleased to hear of it). In particular, consider a Pearson-style system based upon a cubic polynomial. This will be the family of solutions $p(x)$ to the differential equation: $$\frac{dp(x)}{dx}=-\frac{(a+x) p(x)}{c_3 x^3+c_2 x^2+c_1 x+c_0}$$

After some work ... too long to show here, but related to boundary conditions etc ... one ends up with a recurrence relation in raw moments:

$$c_0 (-r) \acute{\mu }_{r-1}-c_2 (r+2) \acute{\mu }_{r+1}-c_3 (r+3) \acute{\mu }_{r+2}-c_1 (r+1) \acute{\mu }_r=-a \acute{\mu }_r-\acute{\mu }_{r+1}$$

Given the boundary condition $\acute{\mu }_0=1$, we enter this recurrence relation into Mathematica as:

enter image description here

Our objective is to find $a$, $c_0$, $c_1$, $c_2$ and $c_3$ in terms of$\acute{\mu }_r$. Putting $r = 0, 1, 2, 3, 4$ yields the required 5 equations (for the 5 unknowns) which we now solve simultaneously to yield the solution:

Z1 = Solve[Table[eqn2[r], {r, 0, 4}], {a, c0, c1, c2, c3}]

The solution is rather long, so I won't print the output here. However, if we work about the mean, taking $\acute{\mu }_1=0$, and $\acute{\mu }_r= \mu_r$ for $r>=2$, the solution reduces to:

enter image description here

which is comparatively compact. The text also provides examples with actual fitting of the third-order (cubic) Pearson, and comparison to the automated second-order (quadratic) Pearson family.

  • $\begingroup$ You set $\mu_0$ to 1, but with $r = 0$, there's a $\mu_{-1}$. Does it have a special value, too? $\endgroup$
    – rcollyer
    Commented Oct 1, 2013 at 13:56
  • $\begingroup$ Hi. The only place $\acute{\mu }_{r-1}$ appears is in the term $c_0 (-r) \acute{\mu }_{r-1}$, and when $r=0$, that whole term is 0 (because of the (-r) multiplier ). $\endgroup$
    – wolfies
    Commented Oct 1, 2013 at 14:03
  • $\begingroup$ You're absolutely correct. I hadn't noticed the $r$ coefficient. Thanks. $\endgroup$
    – rcollyer
    Commented Oct 1, 2013 at 14:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.