I have a matrix with 3 columns. The first columns contains values of a variable $x_1 \in [-1,1]$. The second columns contains values of a variable $x_2 \in [-1,1]$. The third column contains a variable $y$ s.t. $y = p(x_1,x_2)$, where $p$ is the probability density function of a multivariate normal distribution $p = \mathcal{N}(x_1, x_2 | \mu,\Sigma)$, with $\mu \in \mathcal{R}^2$ is the mean of the distribution and $\Sigma \in \mathcal{R}^{2 \times 2}$ is the variance-covariance matrix

How can I estimate $\mu$ and $\Sigma$?

I though that I could estimate $\mu$ by taking the mode of the distribution, that is, $\mu = [\bar{x_1}, \bar{x_2}] = \arg_{\text{max}} p(x_1, x_2)$, this is because the mode of a Gaussian coicides with the mean.

But I have no idea for the variance-covariance matrix.

I have a matrix with 3 columns. The first column contains values of a variable $x_1 \in [-1,1]$. The second column contains values of a variable $x_2 \in [-1,1]$. The third column contains a variable $y$ s.t. $y = p(x_1,x_2)$, where $p$ is the probability density function of a multivariate normal distribution $p = \mathcal{N}(x_1, x_2 | \mu,\Sigma)$, with $\mu \in \mathcal{R}^2$ is the mean of the distribution and $\Sigma \in \mathcal{R}^{2 \times 2}$ is the variance-covariance matrix.

How can I estimate $\mu$ and $\Sigma$?

I thought that I could estimate $\mu$ by taking the mode of the distribution; that is, $\mu = [\bar{x_1}, \bar{x_2}] = \arg_{\text{max}} p(x_1, x_2)$. This is because the mode of a Gaussian coincides with the mean.

But I have no idea for the variance-covariance matrix.