Say I have a random variable $X\sim\mathscr{N}\left(0,\sigma^2\right)$ and another random variable $Y=X+\varepsilon$, where $\varepsilon\sim\mathscr{N}\left(0,\eta^2\right)$ is independent of $X$. I can calculate $Cov\left(X,Y\right)=Cov\left(X,X+\varepsilon\right)=Var(X)=\sigma^2$.
Is there a way to get Mathematica to generate a bivariate normal distribution automatically from the given information? Specifically, rather than working out the covariance matrix myself and then inputting to Mathematica, I would like to be able to do something like (implementing the above example):
\[ScriptCapitalX] = NormalDistribution[0, \[Sigma]]
\[CapitalEpsilon] = NormalDistribution[0, \[Eta]]
\[ScriptCapitalY] = TransformedDistribution[X + \[CurlyEpsilon], {X \[Distributed] \[ScriptCapitalX], \[CurlyEpsilon] \[Distributed] \[CapitalEpsilon]}]
\[Mu] = {Mean[\[ScriptCapitalX]], Mean[\[ScriptCapitalY]]}
\[CapitalSigma] = {{Variance[\[ScriptCapitalX]], Covariance[\[ScriptCapitalX], \[ScriptCapitalY]]}, {Covariance[\[ScriptCapitalX], \[ScriptCapitalY]], Variance[\[ScriptCapitalY]]}}
multi = MultinormalDistribution[\[Mu], \[CapitalSigma]]
want = Expectation[X \[Conditioned] Y = y, {Y, X} \[Distributed] multi]
Obviously, the Covariance function doesn't work that way, but at a minimum I would really like to have Mathematica do the covariance calculation. Also, if there's an easy way to get the mean vector and covariance matrix without explicitly calling Mean, Variance, etc. on individual elements to construct them, that would be even better. I am quite new to the Wolfram language.