Added a test to signify need of transpose of cholesky matrix.

Author: Sayam753

pymc3/tests/test_distributions_random.py

@@ -849,6 +849,12 @@ def ref_rand_chol(size, mu, rowchol, colchol):
                 size, mu, rowcov=np.dot(rowchol, rowchol.T), colcov=np.dot(colchol, colchol.T)
             )
 
+        def ref_rand_chol_transpose(size, mu, rowchol, colchol):
+            colchol = colchol.T
+            return ref_rand(
+                size, mu, rowcov=np.dot(rowchol, rowchol.T), colcov=np.dot(colchol, colchol.T)
+            )
+
         def ref_rand_uchol(size, mu, rowchol, colchol):
             return ref_rand(
                 size, mu, rowcov=np.dot(rowchol.T, rowchol), colcov=np.dot(colchol.T, colchol)
@@ -858,7 +864,7 @@ def ref_rand_uchol(size, mu, rowchol, colchol):
             pymc3_random(
                 pm.MatrixNormal,
                 {"mu": RealMatrix(n, n), "rowcov": PdMatrix(n), "colcov": PdMatrix(n)},
-                size=n,
+                size=100,
                 valuedomain=RealMatrix(n, n),
                 ref_rand=ref_rand,
             )
@@ -867,7 +873,7 @@ def ref_rand_uchol(size, mu, rowchol, colchol):
             pymc3_random(
                 pm.MatrixNormal,
                 {"mu": RealMatrix(n, n), "rowchol": PdMatrixChol(n), "colchol": PdMatrixChol(n)},
-                size=n,
+                size=100,
                 valuedomain=RealMatrix(n, n),
                 ref_rand=ref_rand_chol,
             )
@@ -878,6 +884,22 @@ def ref_rand_uchol(size, mu, rowchol, colchol):
             #     extra_args={'lower': False}
             # )
 
+            # 2 sample test fails because cov becomes different if chol is transposed beforehand.
+            # This implicity means we need transpose of chol after drawing values in
+            # MatrixNormal.random method to match stats.matrix_normal.rvs method
+            with pytest.raises(AssertionError):
+                pymc3_random(
+                    pm.MatrixNormal,
+                    {
+                        "mu": RealMatrix(n, n),
+                        "rowchol": PdMatrixChol(n),
+                        "colchol": PdMatrixChol(n),
+                    },
+                    size=100,
+                    valuedomain=RealMatrix(n, n),
+                    ref_rand=ref_rand_chol_transpose,
+                )
+
     def test_kronecker_normal(self):
         def ref_rand(size, mu, covs, sigma):
             cov = pm.math.kronecker(covs[0], covs[1]).eval()