---
title: "Fitting a continuous non-parametric second-order distribution to data"
canonical: "https://modelassist.epixanalytics.com/space/EA/26575314/Fitting%20a%20continuous%20non-parametric%20second-order%20distribution%20to%20data"
format: markdown
---
When we do not have a great deal of data, a considerable amount of uncertainty will remain about an empirical distribution determined directly from the data. It would be very useful to have the flexibility of using an empirical distribution, i.e. not having to assume a [parametric distribution](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575199/), and also to be able to quantify the uncertainty about that distribution. The following Bayesian technique provides these requirements.
Consider a set of *n* data values {*x**j*} drawn from a distribution, and ranked in ascending order {*x**i*} so *x**i***i**+1*. Data thus ranked are known as the order statistics of {*x*}. We can use these order statistics to construct a second order empirical cumulative distribution
Here are the mathematics behind the technique:
| Individually, each of the values of {*x**j*} may map as a U(0,1) onto the cumulative probability of the parent distribution *F(x)*. We take a U(0,1) distribution as the prior distribution for the cumulative probability for any value of *x*. We can thus use a U(0,1) prior for *P**i* = *F(x**i**)* for the value of the *i*th observation. However, we have the additional information that of *n* values drawn randomly from this distribution, *xi* ranked *i*th, i.e. (*i*-1) of the data values are less than *x**i* and (*n-i*) values are greater than *x**i*. Using [Bayes' Theorem](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26574955) and the [Binomial Theorem](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26574954), the posterior marginal distribution for *P**i* can readily be determined, remembering that *Pi* has a U(0,1) prior and therefore a prior probability density = 1:
> Macro (mathblock)
which is simply the standard Beta distribution Beta(*i*,*n*-*i*+1):
> Macro (mathinline)
(1)
Equation 1 could actually be determined directly from the fact that the beta distribution is the conjugate to the binomial likelihood function and that a U(0,1) = Beta(1,1). The mean of the Beta(*i*, *n-i*+1) distribution equals *i */ (*n*+1): a formula that has been used to estimate the [best-fitting first-order non-parametric cumulative distribution](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575313).
Since *P**i+1*>*P**i*, these Beta distributions are not independent so we need to determine the conditional distribution f(*P**i+1**|P**i*), as follows. The joint distribution f(*P**i**,P**j*) for any two *P**i*, *P**j* is calculated using the Binomial Theorem in a similar manner to the numerator of the equation for f(*P**i**|x**i*; i=1,n), i.e.:
> Macro (mathblock)
where *P**j*>*P**i* and remembering that the prior probability densities for *P**i* and *P**j* equal 1 since they have U(0,1) priors.
For *j*=*i*+1:
> Macro (mathblock)
The conditional probability f(*P**i+1**|P**i*) is thus given by:
> Macro (mathblock)
where *k* is some constant. The corresponding cumulative distribution function F(*P**i+1**|P**i*) is then given by:
> Macro (mathblock)
F(*P**i+1**|P**i*)=1 at *P**i+1*=1, so *k* = (*n*-*i*) and the formula reduces to:
> Macro (mathinline)
(2)
Equations 1 and 2 provide us with the tools to construct a non-parametric second-order distribution for a continuous variable given a data set sampled from that distribution. The distribution for the cumulative probability *P**1* that maps onto the first order statistic X*1* can be obtained from Equation 1 by setting *i* = 1:
> Macro (mathinline)
(3)
The distribution for the cumulative probability *P**2* that maps onto the first order statistic *X**2* can be obtained from Equation 2. F(*P**i+1**|P**i*), being a cumulative distribution function, is Uniform(0, 1) distributed. Thus, writing *U**i+1* to represent a Uniform(0, 1) distribution in place of F(*P**i+1**|P**i*), using the identity 1-U(0, 1) = U(0, 1), and rewriting for *P**i+1*, we obtain:
> Macro (mathinline)
(4)
which gives:
> Macro (mathblock)
> Macro (mathblock)
etc..
Note that each of the *U**2**, U**3**, …U**n* Uniform distributions are independent of each other. |
| --- |
Recapping from the proof, the formulae from Equations 3 and 4 are:
> Macro (mathinline)
(3)
> Macro (mathinline)
(4)
where P*j* represents the estimate of the cumulative distribution function at xj, = F(x*j*). These can be used as inputs to construct a [Cumulative distribution](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575222/), together with subjective estimates of the minimum and maximum values that the variable may take, which can also be assigned subjective distributions.
Model [NonParaCont2](#NonParaCont2) creates a second-order distribution using this technique.
> Macro (anchor)
The links to the NonParaCont2 software specific models are provided here:
Crystal Ball
> Macro (inline-media-image)
[NonParaCont2](/wiki/download/attachments/26579343/NonParaCont2-CB.xlsx?version=1&modificationDate=1506113497000&cacheVersion=1&api=v2)
Note that in Crystal Ball, as F(x) changes from iteration to iteration, the Custom Distribution therefore needs to be capable of [dynamic referencing](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575288/).
##### Using this technique in a second order model
Use [Crystal Ball's 2-D simulation Tool](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575571/), as also illustrated [here](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575534):
- Select "Run" >> "Tools" >> 2D Simulation;
- In Step 1 of 3, select "Forecast" (Cell I14) as target;
- Click Next
- In Step 2 of 3, select the assumption "Distribution (Model)" as Variability (all other assumptions are uncertainty)
- Click Next
- If you run the outer (uncertainty) simulation for 10 trials and run inner (variability) simulation for 1000 trials, this will result in the following second order model:
@Risk
> Macro (inline-media-image)
[NonParaCont2](/wiki/download/attachments/26581760/NonParaCont2-AtRISK.xlsx?version=1&modificationDate=1506097517000&cacheVersion=1&api=v2)
In @Risk the inputs to [RiskCumul](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26581732) can be formulated as: =RiskCumul(min, max, {x*j*}, {P*j*})
##### Using this technique in a second order model
- The uncertainty distributions for *F(x)* are nominated as outputs;
- A smallish number of iterations are run;
- The resultant data are exported back to a spreadsheet;
- Those data are then used to perform multiple simulations (the ”outer loop”) of uncertainty using @RISK’s RiskSimtable function: the ”inner loop” of comes from the Cumulative distribution itself
##### Limitations
There are a few limitations to this technique. In using a [Cumulative distribution](https://epixanalytics.atlassian.net/wiki/spaces/EA/pages/26575222/), one is assuming a histogram style probability distribution function between each of the {x} values. When there are a large number of data points, this approximation becomes irrelevant. However, for small data sets the approximation will tend to accentuate the tails of the distribution: a result of the histogram "squaring-off" effect of using the Cumulative distribution. In other words, the variability will be slightly exaggerated. However, the squaring effect can be reduced, if required, by using some sort of smoothing algorithm and defining points between each observed value. In addition, for small data sets, 'the tails' contribution to the variability will often be more influenced by the subjective estimates of the minimum and maximum values: a fact one can view positively (one is recognizing the real uncertainty about a distribution's tail), and negatively (the smaller the data set, the more the technique relies on subjective assessment).
The fewer the data points, the wider the confidence intervals will become quite naturally and, in general, the more emphasis will be placed on the subjectively defined minimum and maximum values. Conversely, the more data points available, the less influence the minimum and maximum estimates will have on the estimated distribution. In any case, the values of the minimum and maximum only have influence on the width (and therefore height) of the end two histogram bars in the fitted distribution. The fact that the technique is non-parametric, i.e. that no statistical distribution with a particular cumulative distribution function is assumed to be underlying the data, allows the analyst a far greater degree of flexibility and objectivity than that afforded by fitting parametric distributions.
---