Correlation-aware probabilistic data summarization for large-scale multi-block scientific data visualization

Because the maximum likelihood function of the Gaussian copula functions in Eqs. ( 6 ) and ( 7 ) is complex and cannot be solved using partial derivatives, the Broyden–Fletcher–Goldfarb–Shanno (BFGS) method [ 39 ] is used, which is the most popular quasi-Newton algorithm, to accelerate the EM algorithm for this parameter estimation problem with missing data.

We now consider discrete calculation of the copula’s marginal CDF. To calculate the parameters of the Gaussian copula function, it is necessary to pre-calculate the marginal CDFs of value distributions and spatial distributions. The discrete calculation of $H_i$ and $H_i^t(t=1,\mathrm{\cdots },T_n)$ in Eq. ( 6 ) is obvious, i.e., accumulation of probabilities corresponding to the bins satisfying certain conditions. The calculation of $S_i^j$ and $S_i^{j,t}(t=1,\mathrm{\cdots },T_s)$ in Eq. ( 7 ) is also computed discretely. Consider $S_i^j$ as an example: if $\forall p=(x,y,z)\in P_i^j$ , then the value of the corresponding CDF is

(8) $$\begin{array}{c}\hfill S_i^j(p)=\sum {\mathrm{SGmm}}_i^j(p_s)\end{array}$$

where $p_s=(x_s,y_s,z_s)\in P_i^j$ , and $x_s⩽x,y_s⩽y,z_s⩽z$ .

After calculating the CDFs, the parameters of the Gaussian copula functions are estimated. The EM algorithm [ 38 ] is an iterative algorithm used for calculating the maximum likelihood of parameters, and it is widely applied to incomplete data. However, the EM algorithm has a sublinear convergence speed, and the derivative of the $Q$ function in the EM algorithm has no explicit expression. Hence BFGS [ 39 ] is used to accelerate the algorithm. For the case of calculating SC ${}_i^j$ in Eq. ( 7 ), pseudo-code is shown in Algorithm 1.

We determine step size and stopping conditions as follows. The $Q$ function in Algorithm 1 is the expectation of the logarithmic likelihood function of the complete data on the conditional probability distribution of the unobserved data. The step size ${\alpha }_m$ is determined using a backtracking line search so that the energy decreases monotonically. To optimize $g_n({\mathrm{\Theta }}_i^j)$ , ${\alpha }_m$ in the $m\mathrm{th}$ ( $m>1$ ) iteration is initialized to ${\alpha }_{m-1}$ . ${\alpha }_m$ is multiplied by two if the initial ${\alpha }_m$ decreases the energy; ${\alpha }_0=1$ . For the stopping conditions, we set $\eta ={\eta }_1={\eta }_2={10}^{-6}$ and $n_{\max }=m_{\max }=500$ in our experiments.

In Fig. 8 , a schematic diagram of the BFGS-based EM acceleration algorithm for estimating the parameter of a simple bivariate Gaussian copula function is shown to explain the algorithm. The gray curve represents the objective function (the logarithmic likelihood function of ${\mathrm{\Theta }}_i^j$ ). The blue curves represent the lower bound function $g_n({\mathrm{\Theta }}_i^j)$ , which constantly approaches the local optimum, converging to it.

10.1007/s41095-022-0304-6.F008 Fig. 8Schematic diagram of the BFGS-based EM acceleration algorithm.