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.

What I want to do:
1) Define a joint distribution function using CopulaDistribution
2) then transform it with TransformedDistribution
3) to get the new Culmulative Distribution Function (CDF).

The problem:
If I choose correlation 0 or the Product Copula, everything is working perfectly, I obtain a CDF. If I choose a copula different from independence, I don't get a result anymore. In the code below, change the parameter for the copula from 0. to e.g. 0.2, then there's no result for the CDF anymore. (I'm using Mathematica Version 8.0)

The code:

dist = CopulaDistribution[{"Binormal", 0.}, {UniformDistribution[{0, 10}], 
    UniformDistribution[{0, 10}]}]; 

f[x_, y_] = 
  Piecewise[{{x - 3, x > 3 && y > 3}, {x - 6 + y, 
     x > 3 && y <= 3 && x + y - 6 >= 0}}]; 

transdist = 
  TransformedDistribution[f[x, y], {x, y} \[Distributed] dist];

CDF[transdist, 2.]
share|improve this question

2 Answers 2

Copulas (with dependent structures) can result in complicated functional forms. Your example highlights why working with functions as BLACK BOXES can rapidly get one into trouble. It is always a good idea to check both the input and output from such functions at every stage, so that we have some idea of what we are actually working with, and what we could reasonably hope to attain.

For your Binormal copula, for general parameter $a$, the pdf has a complicated form ... check it using PDF[dist, {x, y}] ... which yields the joint pdf of $(X,Y)$ as say $g(x,y)$:

Your transformation equation is:

  f[x_, y_] = Piecewise[{{x-3, x>3 && y>3}, {x-6+y, x>3 && y<=3 && x+y-6 >= 0}}]

You seek the cdf of $Z = f(x,y)$, i.e. $P(f(x,y) < z)$:

where I am using the Prob function from the mathStatica package for Mathematica. Although Mathematica cannot obtain a closed-form, we now have the functional form for the cdf, and it is easy to obtain solutions for given parameter values. For example, here is the solution expressed as a numerical cdf (just swapping Integrate for NIntegrate):

[To zoom in on the above, just open the picture in a new tab or window.]

For your problem, parameter $a=1/5$ and you seek $P(Z<2)$:

 numCDF[2, 1/5]

0.540261

Here is a quick plot of the CDF (when $a=\frac15$):

Plot[numCDF[z, 1/5], {z, -2, 8}, PlotRange -> All, AxesLabel -> {"z", "cdf"}]

It looks almost linear! :) Note the discrete mass at $Z = 0$ that arises due to the structure of your transformation equation f[x,y].

share|improve this answer

I post this just for version 10. It appears to handle this. I would guided by @wolfies as to correctness. Using the definitions from OP:

Plot[Evaluate[CDF[transdist, z]], {z, 0, 10}]

enter image description here

The symbolic definitions:

enter image description here

Note all the following yield 0.545:

CDF[transdist, 2]
Integrate[PDF[transdist, x], {x, 0, 2}] /.HeavisideTheta[x_] -> UnitStep[x]
Probability[f[x, y] < 2, {x, y} \[Distributed] dist]
Probability[z < 2, z \[Distributed] transdist]

and the following yielded: 0.5449

N[Length[Cases[#, _?(# < 2 &)]]/Length[#]] &@RandomVariate[transdist, 100000]

and finally as this just treating as uniform distribution:

reg = ImplicitRegion[
  Reduce[f[x, y] < 2 && 0 < x < 10 && 0 < y < 10, {x, y}], {x, y}];
rp = RegionPlot[reg]
NIntegrate[1, {x, y} \[Element] reg]/100.

The region:

enter image description here

and the integral: 0.545

share|improve this answer

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.