The ML switch

The wCorr package computes correlation coefficients between two vectors of random variables that are jointly bivariate normal. We call the two vectors X and Y.

$$\begin{pmatrix} X \\ Y \end{pmatrix} \sim N \left[ \begin{pmatrix} \mu_x \\ \mu_y \end{pmatrix}, \boldsymbol{\Sigma} \right] $$

where N(A, Σ) is the bivariate normal distribution with mean A and covariance Σ.

Simulation study

To demonstrate the effect of the ML and fast switches a few simulation studies are performed to compare the similarity of the results when the switch is set to TRUE to the result when the switch is set to FALSE. This is done first for the ML switch and then for the fast switch.

Finally, simulations show the implications of these switches on the speed of the computation.

General procedures of the simulation study of unweighted correlations

A simulation is run several times.⁶ For each iteration, the following procedure is used:⁷

• select a true correlation coefficient ρ;
• select the number of observations n;
• generate X and Y to be bivariate normally distributed using a pseudo-Random Number Generator (RNG);
• using a pseudo-RNG, select the number of bins for M and P (t and t′) independently from the set {2, 3, 4, 5};
• select the bin boundaries for M and P (θ and θ′) by sorting the results of (t − 1) and (t′ − 1) draws, respectively, from a normal distribution using a pseudo-RNG;
• confirm that at least 2 levels of each of M and P are occupied (if not, return to the previous step); and
• calculate and record the correlation coefficients.

One of a few possible statistics is then calculated. To compare two levels of a switch the Relative Mean Absolute Deviation is used

$$RMAD= \frac{1}{m} \sum_{j=1}^m | r_{j, \mathtt{TRUE}} - r_{j, \mathtt{FALSE}} | $$

where there are m simulations run, r_j, TRUE and r_j, FALSE are the estimated correlation coefficient for the jth simulated dataset when the switch is set to TRUE and FALSE, respectively. This statistic is called “relative” because it is compared to the other method of computing the statistic, not the true value.

To compare either level to the true correlation coefficient the Root Mean Square Error is used

$$RMSE= \sqrt{ \frac{1}{m} \sum_{j=1}^m (r_j - \rho_j)^2 } $$

where, for the jth simulated dataset, r_j is an estimated correlation coefficient and ρ_j is the value used to generate the data (X, Y, M, and P).

ML switch

A simulation was done using the Cartesian product (all possible combinations of) ML ∈ {TRUE, FALSE}, ρ ∈ (−0.99, −0.95, −0.90, −0.85, ..., 0.95, 0.99), and n ∈ {10, 100, 1000}. Each iteration is run three times to increase the precision of the simulation. The same values of the variables are used in the computation for ML=TRUE as well as for ML=FALSE; and then the statistics are compared between the two sets of results (e.g. ML=TRUE and ML=FALSE).

Figure 1 :. Root Mean Square Error for ML=TRUE and ML=FALSE.

The RMSE for these two options is so similar that the two lines cannot be distinguished for most of the plot. The exact differences are shown for Polychoric in Table 1 : and for Polyserial in Table 2 :. The column labeled, “RMSE difference” shows how much larger the RMSE is for ML=TRUE than ML=FALSE. Because this difference is always positive the RMSE of the ML=FALSE option is always lower. Because of this the ML=TRUE option is preferable only in situations where there is a reason to prefer MLE for some other reason.

Table 1 :. Relative Mean Absolute Deviation between ML=TRUE and ML=FALSE for Polychoric.

Table 2 :. Relative Mean Absolute Deviation between ML=TRUE and ML=FALSE for Polyserial.

For the Polychoric, the agreement between these two methods, in terms of MSE is within 0.001 for n of 10 and decreases to within less than 0.0000 for n of 100 or more. Given the magnitude of these differences the faster method will be preferable.

The final column in the above tables shows the RMAD which compares how similar the ML=TRUE and ML=FALSE results are to each other. Because these values are larger than 0, they indicate that there is not complete agreement between the two sets of estimates. If a user considers the MLE to be the correct estimate then they show the deviation of the ML=FALSE results from the correct results.

fast switch

This section examines the agreement between the pure R implementation of the function that calculates the correlation and the Rcpp and RcppArmadillo implementation, which is expected to be faster. The code can compute with either option by setting fast=FALSE (pure R) or fast=TRUE (Rcpp).

A simulation was done at each level of the Cartesian product of fast ∈ {TRUE, FALSE}, ρ ∈ (−0.99, −0.95, −0.90, −0.85, ..., 0.95, 0.99), and n ∈ {10, 100, 1000}. Each iteration was run 100 times. The same values of the variables are used in the computation for fast=TRUE as well as for fast=FALSE; and then the statistics are compared between the two sets of results.

The plot below shows all differences between the fast=TRUE and fast=FALSE runs for the four types of correlations. Note that differences smaller than 10⁻¹⁶ are indistinguishable from 0 by the machine. Because of this, all values were shown as being at least 10⁻¹⁶ so that they could all be shown on a log scale.

Figure 2 :. Relative Mean Absolute Differences between fast=TRUE and fast=FALSE.

The above shows that differences as a result of the fast argument are never expected to be larger than 10⁻⁹ for any type or correlation type. The Spearman never shows any difference that is different from zero and the Pearson show differences just larger than the smallest observable difference when using double precision floating point values (about 1 × 10⁻¹⁶). This indicates that the computation differences are completely irrelevant for these two types.

For the other two types, it is unclear which one is correct and the agreement that is never more distant than the 10⁻⁹ level indicates that any use that requires precision of less than 10⁻⁹ can use the fast=TRUE argument for faster computation.

Implications for speed

To show the effect of the ML and fast switches on computation a simulation was done at each level of the Cartesian product of ML ∈ {TRUE, FALSE}, fast ∈ {TRUE, FALSE}, ρ ∈ (−0.99, −0.95, −0.90, −0.85, ..., 0.95, 0.99), and n ∈ {10¹, 10^1.25, 10^1.5, ..., 10⁷}. Each iteration is run 80 times when n < 10⁵ and 20 times when n ≥ 10⁵. The same values of the variables are used in the computations at all four combinations of ML and fast. A variety of correlations are chosen so that the results represent an average of possible values of ρ.

The following plot shows the mean computing time (in seconds) versus n.

Figure 3 :. Computation time comparison.

In all cases setting the ML option to FALSE and the fast option to TRUE speeds up–or does not slow down computation. Users wishing for the fastest computation speeds will use ML=FALSE and fast=TRUE.

For the Polychoric, when n is a million observations (n = 10⁷), the speed of a correlation when fast=FALSE is 761 seconds and when the fast=TRUE it is and 8 seconds. When fast=TRUE, setting ML=FALSE speeds computation by 3 seconds.

For the Polyserial, when n is a million observations, the speed of a correlation when ML=TRUE is 652 seconds and when the ML=FALSE it is and 77 seconds. When ML=FALSE, setting fast=TRUE speeds computation by 24 seconds.

#Conclusion Overall the simulations show that the ML option is not more accurate but does add computation burden.

The fast=TRUE and fast=FALSE option are a Rcpp version of the correlation code and an R version, respectively and agree with each other–the differences are not expected to be larger than 10⁻⁹.

Thus users wishing for fastest computation speeds and accurate results can use ML=FALSE and fast=TRUE.