Issues-with-covaraince-matrix.md

@@ -28,10 +28,17 @@ To resolve this issue, we propose to use the Moore-Penrose pseudo-inverse method
 <img src="uploads/4692d03ff868fcc538ab06803102a7c3/pinv.png"  width="440" ><img src="uploads/8d64ddc718f1fe06d54b62ca6028d123/pinv_off_diag.png"  width="440" >
 In the above figures, the quantity 'pinv' refers to the `numpy.linalg.pinv` function. The figures indicate that the inner product using the right-hand inverse matrix is consistent with the left-hand inverse matrix, and also the off-diagonal elements are nearly equal to zero. Therefore, `pinv` more robust.
 
-# Ill-conditioned covariance matrix
+# ill-conditioned covariance matrix
 Now, we discuss that the covariance matrix is ill-conditioned (i.e., nearly singular). Let us focus on the eigenvalues of the covariance matrix.
 <img src="uploads/b4604ec7877c560d6f7aa9232ffd4c08/eigenvalues1.png"  width="440" >
 
 The above plot shows the array of the eigenvalues. This plot indicates that the covariance matrix is ill-conditioned since several eigenvalues are close to zero. Also, the eigenvalues are oscillating around and after the index number 400 (the eigenvalues below ~1e-11), which occurs due to the numerical inaccuracies. This can affect the inverse calculation since the inverse is proportional to the determinant of that matrix. Therefore, we impose a threshold on eigenvalues to calculate the inverse in the `pinv` method. In other words, a reduced singular matrix can be used in the SVD based inverse computation, which can capture the maximum feature of the covariance matrix and also excludes the numerical instabilities. The parameter `rcond` in `numpy.linalg.pinv` function controls this criterion such that singular values less than or equal to `rcond * largest_singular_value` are set to zero.
 
+# Choice of the rank of the covariance matrix
+
+First, we demonstrate an example from the GW190814 event to find an appropriate cut-off on the low-rank approximation of the covariance matrix, and its effect on $`\beta`$. We choose top 5 PE samples and plot $`\beta`$ as a function of the rank of the covariance matrix. The spectrum of the eigenvalues is also plotted. 
+
+<img src="uploads/c93fa55e8c322397cdd710513d862054/beta_lambda_all.png"  width="440" > 
+
+We can see that the median $`\beta`$ (black trace) saturates after a certain number of the top basis vectors have been used to approximate the covariance matrix.  If we use the top 250 basis vectors (rank = 250, instead of the full rank), the Frobenius norm of the residual is $`1.4\times 10^{-8}`$, indicating the extent to which we capture the features in the covariance matrix. The blue trace is for the relative magnitude of the eigenvalues. If we include the eigenvectors with relative importance less than few times 1e-10, you can see that the variance of the $`\beta`$ values starts flaring up again due to errors in the inverse of the covariance matrix. So our suggestion is to use top-p = 250 eigenvectors for the pseudo inverse of the covariance matrix for this event. We feel that this should be a tunable number based on plots like this.