Mean and covariance of a mixture of two distributions This is a homework problem. If I could have some hints on how to solve it that would be great.
Mean and Covariance of a Mixture of Two Distributions: 
$$ I = \left\{ 
  \begin{array}{l l}
    1 & \quad \text{if subject in group 1}\\
    0 & \quad \text{if subject in group 2}
  \end{array} \right.$$
Suppose $\mathbf{Z}=(Z_1,\ldots,Z_p)^T$ is a vector of measurements on a subject with mean and covariances: $$ E[\mathbf{Z}|I=1]=\mathbf{\mu}_1, \quad E[\mathbf{Z}|I=0]=\mathbf{\mu}_2\\Cov[\mathbf{Z}|I=1]=\mathbf{\Sigma}_1, \quad Cov[\mathbf{Z}|I=0]=\mathbf{\Sigma}_2$$
Let $p=P(I=1)$. Find $E[\mathbf{Z}]$ and $Cov[\mathbf{Z}]$. 
My Attempt:
(This makes sense to me) $E[\mathbf{Z}] = E[E[\mathbf{Z}|I]] = E[p\mathbf{\mu}_1+(1-p)\mathbf{\mu}_2] = p\mathbf{\mu}_1+(1-p)\mathbf{\mu}_2$.
(Not sure about this one) $$ \begin{align*}Cov[\mathbf{Z}] &= Cov[E[\mathbf{Z}|I]] + E[Cov[\mathbf{Z}|I]]\\&= Cov[p\mathbf{\mu}_1+(1-p)\mathbf{\mu}_2] + E[Cov[\mathbf{Z}|I]]\\&= ???\end{align*}$$
 A: Here is my solution if anyone is curious (sorry it is so ugly):
We begin by finding the expected value using iterated expectation
$\begin{align*} E[\mathbf{Z}]&=E[E[\mathbf{Z}|I]]  \\ &= E[p\mu_1+(1-p)\mu_2]\\ &= p\mu_1+(1-p)\mu_2  \end{align*}$
Next we note that the group specific covariances can be written as
$\begin{align*} Cov[\mathbf{Z}|I=1] &= E[\mathbf{Z}\mathbf{Z}^T|I=1]-\mu_1\mu_1^T = \mathbf{\Sigma}_1\\   Cov[\mathbf{Z}|I=0] &= E[\mathbf{Z}\mathbf{Z}^T|I=0]-\mu_2\mu_2^T = \mathbf{\Sigma}_2\\   E[\mathbf{Z}\mathbf{Z}^T|I=1] &= \mathbf{\Sigma}_1 + \mu_1\mu_1^T\\   E[\mathbf{Z}\mathbf{Z}^T|I=0] &= \mathbf{\Sigma}_2 + \mu_2\mu_2^T  \end{align*} $
Using this same equation for covariance we get
$\begin{align*}  Cov[\mathbf{Z}] &= E[\mathbf{Z}\mathbf{Z}^T] - E[\mathbf{Z}]E[\mathbf{Z}]^T\\  &=E[E[\mathbf{Z}\mathbf{Z}^T|I]] - (p\mu_1+(1-p)\mu_2)(p\mu_1+(1-p)\mu_2)^T\\  &= E[p(\mathbf{\Sigma}_1 + \mu_1\mu_1^T)+(1-p)(\mathbf{\Sigma}_2 + \mu_2\mu_2^T)] - (p\mu_1+(1-p)\mu_2)(p\mu_1^T+(1-p)\mu_2^T)\\  &= p\mathbf{\Sigma}_1 + p\mu_1\mu_1^T+(1-p)\mathbf{\Sigma}_2 + (1-p)\mu_2\mu_2^T - (p\mu_1+(1-p)\mu_2)(p\mu_1^T+(1-p)\mu_2^T)\\  &=p\mathbf{\Sigma}_1 + p\mu_1\mu_1^T+(1-p)\mathbf{\Sigma}_2 + (1-p)\mu_2\mu_2^T - (p^2\mu_1\mu_1^T+p(1-p)\mu_1^T\mu_2+p(1-p)\mu_1\mu_2^T+(1-p)^2\mu_2\mu_2^T)\\  &=p\mathbf{\Sigma}_1+(1-p)\mathbf{\Sigma}_2+(p-p^2)(\mu_1\mu_1^T)-p(1-p)\mu_1^T\mu_2-p(1-p)\mu_1\mu_2^T+(p-p^2)(\mu_2\mu_2^T)\\  &=p\mathbf{\Sigma}_1+(1-p)\mathbf{\Sigma}_2+p(1-p)\left[(\mu_1\mu_1^T)-\mu_1^T\mu_2-\mu_1\mu_2^T+(\mu_2\mu_2^T)\right]\\  &= p\mathbf{\Sigma}_1+(1-p)\mathbf{\Sigma}_2+p(1-p)(\mu_1-\mu_2)(\mu_1^T-\mu_2^T)\\  &= p\mathbf{\Sigma}_1+(1-p)\mathbf{\Sigma}_2+p(1-p)(\mu_1-\mu_2)(\mu_1-\mu_2)^T\\  \end{align*}$
