fix(riskmodel): keep structured covariance on one sample denominator

examples/portfolio/README.md

@@ -33,6 +33,7 @@ We use China stock market data for our example.
    ```
 
 Here we use a **Statistical Risk Model** implemented in `qlib.model.riskmodel`.
+`StructuredCovEstimator` writes factor exposures, a factor covariance, and specific risk (the square root of residual variance). The factor covariance and the residual variance both divide by `N - 1`, the same convention as `numpy.cov`, and a one-factor model still writes a 1×1 factor covariance. For a PCA fit the diagonal matches the sample variances of the window. Factor analysis does not have that identity.
 However users are strongly recommended to use other risk models for better quality:
 * **Fundamental Risk Model** like MSCI BARRA
 * [Deep Risk Model](https://arxiv.org/abs/2107.05201)

qlib/model/riskmodel/structured.py

@@ -21,6 +21,14 @@ class StructuredCovEstimator(RiskModel):
     Therefore, the structured covariance can be estimated by
         cov(X.T) = F @ cov(B.T) @ F.T + diag(var(U))
 
+    Both moments divide by `N - 1`, the same convention as `numpy.cov`.
+    That common denominator is not an unbiased specific variance: a k-factor
+    fit has already used k degrees of freedom. For a PCA fit the factors are
+    orthogonal to the residuals, so the diagonal of this estimate equals the
+    sample variances of `X`. Factor analysis does not have that identity.
+    `cov(B.T)` is always `(num_factors, num_factors)`. A single factor is
+    included: `numpy.cov` would otherwise return a scalar.
+
     In finance domain, there are mainly three methods to design `F` [1][2]:
         - Statistical Risk Model (SRM): latent factor models major components
         - Fundamental Risk Model (FRM): human designed factors
@@ -83,8 +91,11 @@ def _predict(self, X: np.ndarray, return_decomposed_components=False) -> Union[n
         F = model.components_.T  # variables x factors
         B = model.transform(X)  # observations x factors
         U = X - B @ F.T
-        cov_b = np.cov(B.T)  # factors x factors
-        var_u = np.var(U, axis=0)  # diagonal
+        # One factor is a 1-d sample, and np.cov then returns a scalar.
+        cov_b = np.atleast_2d(np.cov(B.T))  # factors x factors, sample ddof
+        # np.var defaults to dividing by N. That is smaller than np.cov's N - 1
+        # denominator, so specific variance would not match the factor block.
+        var_u = np.var(U, axis=0, ddof=1)  # diagonal
 
         if return_decomposed_components:
             return F, cov_b, var_u

tests/test_structured_cov_estimator.py

@@ -106,6 +106,42 @@ def test_decomposition(self):
 
         self.assertTrue(if_identical)
 
+    def test_single_factor_covariance_is_a_matrix(self):
+        # np.cov of one factor is a scalar. The estimator must still return a
+        # square covariance and a 1x1 factor-covariance block.
+        rng = np.random.RandomState(0)
+        n_obs, n_var = 40, 6
+        X = rng.normal(size=(n_obs, n_var))
+
+        for factor_model in ("pca", "fa"):
+            estimator = StructuredCovEstimator(
+                factor_model=factor_model,
+                num_factors=1,
+                scale_return=False,
+                assume_centered=True,
+            )
+            cov = estimator.predict(X, is_price=False)
+            self.assertEqual(cov.shape, (n_var, n_var))
+            self.assertTrue(np.allclose(cov, cov.T))
+
+            F, cov_b, var_u = estimator.predict(X, is_price=False, return_decomposed_components=True)
+            self.assertEqual(F.shape, (n_var, 1))
+            self.assertEqual(cov_b.shape, (1, 1))
+            self.assertEqual(var_u.shape, (n_var,))
+
+    def test_pca_diagonal_matches_sample_variance(self):
+        # PCA scores are orthogonal to the residuals, so the structured
+        # covariance must reproduce sample variances on its diagonal.
+        rng = np.random.RandomState(1)
+        n_obs, n_var, n_factor = 25, 5, 2
+        X = rng.normal(size=(n_obs, n_var))
+        estimator = StructuredCovEstimator(num_factors=n_factor, scale_return=False, assume_centered=False)
+
+        cov = estimator.predict(X, is_price=False)
+        sample = np.cov(X.T)
+
+        np.testing.assert_allclose(np.diag(cov), np.diag(sample), atol=1e-8)
+
 
 if __name__ == "__main__":
     unittest.main()