$\newcommand{\tr}{\operatorname{tr}}$For a matrix-variate Gaussian distribution, the negative log marginal likelihood is
\begin{equation}\label{matrixLikelihood} \mathcal{L} = \frac{nd}{2}\ln(2\pi) + \frac{d}{2}\ln \det(K') + \frac{n}{2}\ln \det(\Omega) + \frac{1}{2}\tr((K')^{-1}Y\Omega^{-1}Y^{\mathrm{T}}), \end{equation} where $n$ is the number of input points, d is the dimension of output $Y \in \mathbb{R}^{n \times d}$, $K' = K + \sigma^2_n I$, $K \in \mathbb{R}^{n \times n }$ is the column covariance matrix of matrix-variate distribution with many kernel parameters $\theta_i$, for example, squared exponential (SE) kernel, $K = [k_{ij}], k_{ij} = \theta_1 \exp( - \frac{- (x_i - x_j)^2}{2 \theta_2^2})$, $\Omega \in \mathbb{R}^{d \times d}$ is row covariance matrix. In addition, both $K$ and $\Omega$ are postive semi-definite, and thus we can write $\Omega = \Phi \Phi^{\mathrm{T}}$, where $$ \Phi = \left[ \begin{matrix} \phi_{11} & 0 & \cdots & 0 \\ \phi_{21} & \phi_{22} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ \phi_{d1} & \phi_{d2} & \cdots & \phi_{dd} \\ \end{matrix} \right]. $$ To guarantee the uniqueness of $\Phi$, the diagonal elements are restricted to be positive and denote $\varphi_{ii} = \ln(\phi_{ii})$ for $i = 1,2,\cdots,d$.
How to prove the existence of a maximum likelihood estimator (MLE) with respect to parameter $\sigma_n^2$, $\theta_i$, $\phi_{ij}$ and $\varphi_{ii}$. Here our parameter space are $\mathbb{R}^n_{+}$ for $\sigma_n^2$, $\theta_i$, and $\varphi_{ii}$ and $\mathbb{R}^n$ for $\phi_{ij}$.
I was trying to use the theorem that the continuous function on a compact set must have minimum and maximum. But the problem is now our parameter space is not compact. Therefore, I tried to check the monotonicity of the likelihood function at infinity, so I obtain the derivative: with respect to parameter $\sigma_n^2$, $\theta_i$, $\phi_{ij}$ and $\varphi_{ii}$ are as follows \begin{align*} \frac{\partial \mathcal{L}}{\partial \sigma_n^2} &= \frac{d}{2}\tr((K')^{-1}) - \frac{1}{2}\tr(\alpha_{K'}\Omega^{-1}\alpha_{K'}^{\mathrm{T}}),\\ \frac{\partial \mathcal{L}}{\partial \theta_i} & = \frac{d}{2}\tr\left((K')^{-1}\frac{\partial K_{\theta}}{\partial \theta_i}\right) - \frac{1}{2}\tr\left(\alpha_{K'}\Omega^{-1}\alpha_{K'}^{\mathrm{T}}\frac{\partial K_{\theta}}{\partial \theta_i}\right), \\ \frac{\partial \mathcal{L}}{\partial \phi_{ij}} & = \frac{n}{2}\tr[\Omega^{-1}(\mathbf{E}_{ij}\Phi^{\mathrm{T}} + \Phi \mathbf{E}_{ij})] - \frac{1}{2}\tr[\alpha_{\Omega}(K')^{-1}\alpha_{\Omega}^{\mathrm{T}}(\mathbf{E}_{ij}\Phi^{\mathrm{T}} + \Phi \mathbf{E}_{ij})], \\ \frac{\partial \mathcal{L}}{\partial \varphi_{ii}} & = \frac{n}{2}\tr[\Omega^{-1}(\mathbf{J}_{ii}\Phi^{\mathrm{T}} + \Phi \mathbf{J}_{ii})] - \frac{1}{2}\tr[\alpha_{\Omega}(K')^{-1}\alpha_{\Omega}^{\mathrm{T}}(\mathbf{J}_{ii}\Phi^{\mathrm{T}} + \Phi \mathbf{J}_{ii})], \end{align*} where $\alpha_{K'} = (K')^{-1}Y$, $\alpha_{\Omega} = \Omega^{-1}Y^{\mathrm{T}}$, $\mathbf{E}_{ij}$ is the $d \times d$ elementary matrix having unity in the (i,j)-th element and zeros elsewhere, and $\mathbf{J}_{ii}$ is the same as $\mathbf{E}_{ij}$ but with the unity being replaced by $e^{\varphi_{ii}}$.
But I don't know how to evaluate these derivatives positive or negative.
Similarly, is it possible to prove the MLE existence for matrix-variate student-t distribution \begin{eqnarray*}\label{MultiLikelihoodT} % \nonumber % Remove numbering (before each equation) \mathcal{L} &=& \frac{1}{2}(\nu+ d+n -1) \ln \det(\mathbf{I}_n + (K')^{-1}Y\Omega^{-1}Y^{\mathrm{T}}) + \frac{d}{2}\ln \det(K') + \frac{n}{2}\ln \det(\Omega) \\ & & {} + \ln\Gamma_n \left(\frac{1}{2}(\nu + n -1)\right) - \ln \Gamma_n \left(\frac{1}{2}(\nu + d + n -1)\right) + \frac{1}{2}dn\ln\pi \\ &=& \frac{1}{2}(\nu+ d+n -1) \ln \det(K' +Y\Omega^{-1}Y^{\mathrm{T}}) - \frac{\nu + n -1}{2}\ln \det(K') \\ & & {} + \ln\Gamma_n \left(\frac{1}{2}(\nu + n -1)\right) - \ln \Gamma_n \left(\frac{1}{2}(\nu + d + n -1)\right)+ \frac{n}{2}\ln \det(\Omega)+ \frac{1}{2}dn\ln\pi, \end{eqnarray*} where $\nu$ is the degree of freedom of student-t distribution. \begin{align*} \frac{\partial \mathcal{L}}{\partial \nu} &= \frac{1}{2}\ln \det(U) - \frac{1}{2}\ln \det(K') + \frac{1}{2}\psi_n(\frac{1}{2}\tau) - \frac{1}{2}\psi_n\left(\frac{1}{2}(\tau + d)\right), \\ \frac{\partial \mathcal{L}}{\partial \sigma^2_n} & = \frac{(\tau+d)}{2}\tr(U^{-1}) - \frac{\tau}{2}\tr((K')^{-1}), \\ \frac{\partial \mathcal{L}}{\partial \theta_i} & = \frac{(\tau+d)}{2}\tr\left(U^{-1} \frac{\partial K_{\theta}}{\partial \theta_i}\right) - \frac{\tau}{2}\tr\left(\Sigma^{-1} \frac{\partial K_{\theta}}{\partial \theta_i}\right),\\ \frac{\partial \mathcal{L}}{\partial \phi_{ij}} & = - \frac{(\tau +d)}{2}\tr[U^{-1}\alpha_{\Omega}^{\mathrm{T}}(\mathbf{E}_{ij}\Phi^{\mathrm{T}} + \Phi \mathbf{E}_{ij})\alpha_{\Omega}] + \frac{n}{2}\tr[\Omega^{-1}(\mathbf{E}_{ij}\Phi^{\mathrm{T}} + \Phi \mathbf{E}_{ij})], \\ \frac{\partial \mathcal{L}}{\partial \varphi_{ii}} & = -\frac{(\tau +d)}{2}\tr[U^{-1}\alpha_{\Omega}^{\mathrm{T}}(\mathbf{J}_{ii}\Phi^{\mathrm{T}} + \Phi \mathbf{J}_{ii})\alpha_{\Omega}] + \frac{n}{2}\tr[\Omega^{-1}(\mathbf{J}_{ii}\Phi^{\mathrm{T}} + \Phi \mathbf{J}_{ii})], \end{align*} where $U = K' + Y\Omega^{-1}Y^{\mathrm{T}}$, $\tau = \nu + n -1$ and $\psi_n(\cdot)$ is the derivative of the function $\ln \Gamma_n(\cdot)$ with respect to $\nu$.