Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I am using some extreme value fitting method which results in a parametric distribution for values exceeding some threshold, all values $\geq 0$.

For smaller values I'd like to use a smooth kernel distribution.

Pasting these two together, there is usually a small gap or jump in the distribution. So I looked for the possibility to 'fix' a point in the smoothed distribution, but found no such thing. Is there any way 'force' the smooth kernel distribution density to take a right end value?

It doesn't make a big difference, but there is simply no real reason for such a jump in the density.

Update: Of course it is possible to shift the empirical smoothed density and then scale the empirical density part to obtain a continuous density, but I am not really convinced. To me, it does not seem like the optimal solution. I'll attach two plots to illustrate the problem. enter image description here enter image description here Thank you!

share|improve this question
    
Just curious how you use the SmoothKernelDistribution on your data and how you "paste them together". If you do this incorrectly you'll end up with a distribution that doesn't add up to a probability of 1. Could you provide a small example? –  Sjoerd C. de Vries Nov 21 '12 at 15:24
    
Maybe relevant: there's a bug in SmoothKernelDistribution that I reported, see the discussion after @rm -rf's answer to How to create a heatmap from list of coordinates? –  Jens Nov 21 '12 at 15:30
    
@SjoerdC.deVries: The probability of piercing the threshold is some $0<c<1$. I simply scale the empirical (smoothed) distribution so it adds up to a probability measure. –  user13655 Nov 21 '12 at 16:04
add comment

1 Answer

up vote 3 down vote accepted

Here is a possible solution.

Given a cut point c we can create a kernel density estimate kde truncated on {-Infinity, c} and a tail estimate tail truncated on {c, Infinity} with the additional restriction that the PDF of the tail at c must be equal to PDF[kde, c].

To solve for this point in a way that is robust enough for a variety of data sets and distributions we need to set up special rules for dealing with likelihoods.

The function bound will be used to ensure that indeterminate values of the LogLikelihood for our restricted tail will allow optimization functions to keep searching.

bound[Infinity] := $MinMachineNumber 
bound[-Infinity] := $MinMachineNumber
bound[n_?NumericQ] := n

Now for the tail estimator. The 1000 $MachineEpsilon is to nudge us away from restricted boundaries (and may need adjusted). The initial call to FindDistributionParameters attempts to get good starting values.

tailEstimate[dist_, data_, prob_, c_] := 
 Block[{params, est, tdist},
  params = FindDistributionParameters[data, dist];
  tdist = TruncatedDistribution[{c, \[Infinity]}, dist];
  est = Apply[
    FindMaximum, {{bound[
       LogLikelihood[tdist, Cases[data, x_ /; x > c]]], 
      DistributionParameterAssumptions[tdist] && 
       PDF[tdist, c + 1000 $MachineEpsilon] == prob}, 
     params /. Rule -> List}];
  tdist /. est[[2]]
  ]

Now to put it together with the kernel density estimator.

kde[data_, bw_, c_, tailModel_] :=
 Block[{kern, tail},
  kern = TruncatedDistribution[{-\[Infinity], c}, KernelMixtureDistribution[data, bw]];
  tail = tailEstimate[tailModel, data, PDF[kern, c], c];
  MixtureDistribution[{1, 1}, {kern, tail}]
  ]

Here is an example using heavy tailed data.

data = RandomVariate[ParetoDistribution[1, .5], 1000];
mdist = kde[data, .5, 5, ParetoDistribution[a, b]];
Plot[PDF[mdist, t], {t, 0, 20}, PlotRange -> {0, .2}]

enter image description here

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.