I know how to prove covariance by the linearity of expectation end up with $E(XY)-E(X)E(Y)$, but I want to separate the case to digest deeply, I am a bit stuck at a discrete case below. How should I play around with this joint pmf.
$cov(X,Y)=E[(X-\mu_X)(Y-\mu_Y)]=\sum_{i=1}^k\sum_{j=1}^l(x_j-\mu_X)(y_i-\mu_Y)P(X=x_j,Y=y_i)$
Much appreciated for discrete and continuous case separately.
I don't think that I understand what your entire question, but I do understand that you want to know how to proceed from $\sum_{i=1}^k\sum_{j=1}^l(x_j-\mu_X)(y_i-\mu_Y)P(X=x_j,Y=y_i)$ for a direct proof of the covariance formula $\mathrm{cov}(X,Y) = E(XY) - E(X)E(Y)$. To that end, note that \begin{align} \sum_{i=1}^k\sum_{j=1}^l(x_j-\mu_X)&(y_i-\mu_Y)P(X=x_j,Y=y_i) \\&= \sum_{i=1}^k\sum_{j=1}^l[x_j y_i - \mu_X y_i - \mu_Y x_j + \mu_X \mu_Y]P(X=x_j,Y=y_i) \\ & = \sum_{i=1}^k\sum_{j=1}^lx_j y_iP(X=x_j,Y=y_i) - \mu_X\sum_{i=1}^k\sum_{j=1}^l y_iP(X=x_j,Y=y_i) \\ & \quad - \mu_Y \sum_{i=1}^k\sum_{j=1}^lx_jP(X=x_j,Y=y_i) \\ & \quad +\mu_X\mu_Y \sum_{i=1}^k\sum_{j=1}^l P(X=x_j,Y=y_i) \\ & = \sum_{i=1}^k\sum_{j=1}^lx_j y_iP(X=x_j,Y=y_i) - \mu_X\sum_{i=1}^k y_iP(Y=y_i) \\ & \quad - \mu_Y \sum_{j=1}^lx_jP(X=x_j) +\mu_X\mu_Y \\ & = E(XY) - \mu_X E(Y) - \mu_YE(X) + \mu_X\mu_Y \\ & = E(XY) - E(X)E(Y) - E(Y)E(X) + E(X)E(Y) \\ & = E(XY) - E(X)E(Y). \end{align} The proof for the continuous is essentially the same except that the probability mass function is replaced by a probability density function and the iterated sums are replaced with iterated integrals.