Calculates a correlation matrix \(\mathbf{R}\) or covariance matrix \(\mathbf{\Sigma}\) using the following steps:
1) An input matrix \(\mathbf{X}\) with \(N\) columns is split into \(k\) equal size blocks (+ a possible remainder block) \(A_1, A_2, \ldots, A_k\) of size \(n\). The block size can be defined by the user, size = 2000 is a good value because cor can handle this quite quickly (~ 400 ms). For example, if the matrix has 13796 columns, the split will be \(A_1 = 1 \ldots 2000; A_2 = 2001 \ldots 4000; A_3 = 4001 \ldots 6000; A_4 = 6000 \ldots 8000 ; A_5 = 8001 \ldots 10000; A_6 = 10001 \ldots 12000; A_7 = 12001 \ldots 13796\).
2) For all pairwise combinations of blocks \(k \choose 2\), the \(n \times n\) correlation sub-matrix is calculated. If y = NULL, \(\mathrm{cor}(A_1, A_1), \mathrm{cor}(A_1, A_2), \ldots, \mathrm{cor}(A_k, A_k)\), otherwise \(\mathrm{cor}(A_1, y), \mathrm{cor}(A_2, y), \ldots, \mathrm{cor}(A_k, y)\).
3) The sub-matrices are transferred into a preallocated \(N \times N\) empty matrix at the corresponding position (where the correlations would usually reside). To ensure symmetry around the diagonal, this is done twice in the upper and lower triangle. If y was supplied, a \(N \times M\) matrix is filled, with \(M\) = number of columns in y.
Since the resulting matrix is in 'hdf5r' format, one has to subset to extract regions into normal matrix-like objects. See 'Examples'.