1
$\begingroup$

Consider a random vector {s,c} with a bivariate normal distribution. For a vector of positive scalars {a, ß, σz}, I'm interested in calculating (numerically) the probability

NProbability[c < (1 - CDF[NormalDistribution[ ß*(s - c), σz], a]),{s,c} \[Distributed] BinormalDistribution[{μs, μc}, {σs, σc}, ρ]]

Is there a way to write this same calculation using only NIntegrate?

What I've done so far

I've tried re-writing the probability, solving for s on one side of the inequality, and nesting the integrals:

f1[c_?NumericQ,μs_, μc_, σs_, σc_, σz_, ρ_, a_, ß_]:=NIntegrate[PDF[BinormalDistribution[{μs, μc}, {σs, σc}, ρ],{s,c}],{s, c + (a-σz*InverseCDF[NormalDistribution[],1-c])(ß)^-1,\[Infinity]},]


f2[μs_, μc_, σs_, σc_, σz_, ρ_, a_, ß_]:=NIntegrate[f1[c,μs, μc, σs, σc, σz, ρ, a, ß],{c,-\[Infinity],\[Infinity]}]

This approach unfortunately doesn't work because the computation gets stuck with InverseCDF[NormalDistribution[],1-c] for c below zero or above one.

Parameter values

The scalars and distribution parameters are not important. Here is a starting set of values that can be used for reference:

{μs, μc, σs, σc, σz, ρ, a, ß} = {.35, .5, 1.1, 1.2, 1.3, .25, 1, .5}

Improve this question
$\endgroup$
5
  • 2
    $\begingroup$ Do you have values for σs, σc, ρ, ß, a, σz, μs, μc... etc? NIntegrate expects these to be numerical - as does NProbability. $\endgroup$
    – flinty
    Jun 22, 2020 at 23:48
  • $\begingroup$ @flinty Yes, I do have values for these. I get NProbability to give me the required calculation - but I’d like to reproduce that calculation with NIntegrate $\endgroup$
    – OO_SE
    Jun 23, 2020 at 0:11
  • 1
    $\begingroup$ Can you give the values then? The code in your NIntegrate doesn't work or has been formatted wrongly. I get ... not a valid limit of integration $\endgroup$
    – flinty
    Jun 23, 2020 at 0:42
  • $\begingroup$ Thanks, I've edited the question and added an example set of parameter values, in case that's helpful. The code for NIntegrate, as written above, doesn't work - as you say, it rejects the limit of integration. But if you write this up nesting two integrals (one over s and a second one over c) then you see that if fails due to c<0 and c>1 $\endgroup$
    – OO_SE
    Jun 23, 2020 at 1:11
  • $\begingroup$ I'll expand the question, spelling out the nested integrals, to make my point clear $\endgroup$
    – OO_SE
    Jun 23, 2020 at 1:14

1 Answer 1

Reset to default
1
$\begingroup$ 
{μs, μc, σs, σc, σz, ρ, a, ß} = 
   {.35, .5, 1.1, 1.2, 1.3, .25, 1, .5} // Rationalize;

ineq = c < (1 - CDF[NormalDistribution[ß*(s - c), σz], a]) // 
  FullSimplify

(* 2 c < Erfc[(5 (2 + c - s))/(13 Sqrt[2])] *)

NIntegrate is an available Method for NProbability

NProbability[ineq,
 {s, c} \[Distributed] 
  BinormalDistribution[{μs, μc}, {σs, σc}, ρ],
 WorkingPrecision -> 30,
 Method -> {"NIntegrate",
   {MinRecursion -> 15, MaxRecursion -> 25}}]

enter image description here

(* 0.419500831140737615758538073412 *)

However, the "MonteCarlo" method does not result in warning messages.

NProbability[ineq,
 {s, c} \[Distributed] 
  BinormalDistribution[{μs, μc}, {σs, σc}, ρ],
 WorkingPrecision -> 30, 
 Method -> {"MonteCarlo", "SamplingIncrement" -> 10^4}]

(* 0.415478873239436619718309859154930 *(
Improve this answer
$\endgroup$
1
  • $\begingroup$ Thanks for this. I do get NProbability to give me the result, as you point out her, but what I'd like is to be able to reproduce this writing down the calculation explicitly in terms of NIntegrate, rather than the NProbability function $\endgroup$
    – OO_SE
    Jun 23, 2020 at 2:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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