2
$\begingroup$

Suppose I have a bivariate normal random vector $(Y,Z)'$ where $Y\sim N(μ1,\sigma_1^2)$, $Z\sim N(μ2,\sigma_2^2)$, and $Cov[Y,Z]=\sigma_{12}$. If I want Mathematica to compute the conditional expectation of $Y|Z=z$, then the following command

Expectation[Y \[Conditioned] Z == z, {Y, Z} \[Distributed] MultinormalDistribution[{μ1, μ2}, {{Subscript[\[Sigma],1]^2, Subscript[\[Sigma], 12]}, {Subscript[\[Sigma], 12], Subscript[\[Sigma], 2]^2}}]]

returns the correct answer:

$$E[Y|Z] = \mu1+ \frac{ \sigma_{12} }{\sigma_2^2}(z- \mu2).$$

Is there a way to have Mathematica compute $E[Y| |Z|=z]$? I have tried

Expectation[Y \[Conditioned] Abs[Z] == z, {Y, Z} \[Distributed] MultinormalDistribution[{μ1, μ2}, {{Subscript[\[Sigma],1]^2, Subscript[\[Sigma], 12]}, {Subscript[\[Sigma], 12], Subscript[\[Sigma], 2]^2}}]]

But as output Mathematica returns what I have written as input. More generally, I would like to know whether it is possible to have Mathematica compute $E[Y| |Z|=z, X=x, W=w]$, with $(Y,Z,X,W)$ multivariate normally distributed random variables.

$\endgroup$
3
  • 1
    $\begingroup$ @bills I'm certain that's not what the OP wants. What that gets is just the marginal mean for $Y$. Rather than use Abs[z]==z one could use a==Cos[a] and get the same answer. $\endgroup$
    – JimB
    Nov 13, 2019 at 18:39
  • 2
    $\begingroup$ Your question is not well-defined, and has no solution as stated. It is not sufficient to know that $Y$ is Normal and that $Z$ is Normal and that their covariance is known. You need to know the joint pdf of $(Y,Z)$: that is not stated in the question, and absent which the rest does not follow. $\endgroup$
    – wolfies
    Nov 13, 2019 at 22:39
  • $\begingroup$ @wolfies, you are right. I will amend the question. $\endgroup$
    – Benno
    Nov 14, 2019 at 5:40

1 Answer 1

4
$\begingroup$

If (as stated by @wolfies in a comment) the joint distribution is a bivariate normal, then here's a hand-waving approach (which I think could be made more concrete).

Because you can get the conditional mean for $Z=z$, then the conditional mean for $|Z|=z$ should be a weighted average of the conditional means for $Z=z$ and $Z=-z$ with the weights being the marginal densities of $Z$ at $z$ and $-z$:

dist = MultinormalDistribution[{μ1, μ2}, {{Subscript[σ, 1]^2, Subscript[σ, 12]}, {Subscript[σ, 12], Subscript[σ, 2]^2}}];

μNegative = Expectation[Y \[Conditioned] Z == -z, {Y, Z} \[Distributed] dist] // FullSimplify;
μPositive = Expectation[Y \[Conditioned] Z == z, {Y, Z} \[Distributed] dist] //    FullSimplify;
pdfNegative = PDF[NormalDistribution[μ2, Subscript[σ, 2]], -z];
pdfPositive = PDF[NormalDistribution[μ2, Subscript[σ, 2]], z];
μ = (μNegative pdfNegative + μPositive pdfPositive)/(pdfNegative + pdfPositive) // FullSimplify

Conditional mean

As a somewhat feeble check, consider taking a large random sample from the bivariate distribution and then keep only those values where the realized values of the absolute value of $Z$ are "close" to $z$:

parms = {Subscript[σ, 12] -> 2/3, μ1 -> 2, μ2 -> 2, 
   z -> 1, Subscript[σ, 1] -> 1, Subscript[σ, 2] -> 1};
SeedRandom[12345];
yz = RandomVariate[dist /. parms, 1000000];
yz = Select[yz, (Abs[Abs[#[[2]]] - Abs[z]] /. parms) < 0.01 &];
Length[yz]
(* 5013 *)
Mean[yz[[All, 1]]]
(* 1.31002 *)
μ /. parms // N 
(* 1.30935 *)

Here is a generalization to find

$$E(X_1 |\,\,\, |X_2|=x_2, X_3=x_3, X_4=x_4, \ldots, X_n=x_n)$$

conditionalMean[n_] := Module[{mean, pdf, cm},
  mean = Expectation[X[1] \[Conditioned] (Table[X[i] == x[i], {i, 2, n}] /. List -> And),
    Table[X[i], {i, n}] \[Distributed] MultinormalDistribution[Table[μ[i], {i, n}], 
      Table[σ[Min[i, j], Max[i, j]], {i, n}, {j, n}]]];
  pdf = PDF[MultinormalDistribution[Table[μ[i], {i, 2, n}], 
     Table[σ[Min[i, j], Max[i, j]], {i, 2, n}, {j, 2, n}]], Table[x[i], {i, 2, n}]];
  cm = (mean*pdf + (mean /. x[2] -> -x[2])*(pdf /. x[2] -> -x[2]))/(pdf + (pdf /. x[2] -> -x[2]));
  cm /. σ[i_, i_] -> Subscript[σ, i]^2 /. σ[i_, j_] -> Subscript[σ, i, j] /.
    μ[i_] -> Subscript[μ, i] /. x[i_] -> Subscript[x, i]]

conditionalMean[2] // FullSimplify

Conditional mean with n=2

conditionalMean[3] // FullSimplify

Conditional mean with n=3

So there looks to be a pattern here. conditionalMean[4] shows up quickly but applying FullSimplify afterwards might take forever.

$\endgroup$
2
  • $\begingroup$ Thanks! Neat and really useful! $\endgroup$
    – Benno
    Nov 15, 2019 at 9:26
  • $\begingroup$ Yes, I have tried conditionalMean[4] with FullSimplify and aborted after about one hour $\endgroup$
    – Benno
    Nov 15, 2019 at 11:15

Your Answer

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

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