TransWikia.com

Failure in MultinormalDistribution

Mathematica Asked on December 27, 2020

My RandomVariate fails to evaluate when sampling from a multinormal distribution, any idea what’s happening here?

This is Mathematica 12.1 on Mac:

d = 6;
sigma = {{0.2456589529728737`, 0.06467112508291382`, 
    0.07381608109139866`, 0.01761648518373636`, 0.05925859721052876`, 
    0.045898960400932734`}, {0.06467112508291387`, 
    0.35501955566599797`, 
    0.05242856376555134`, -0.009970312195220573`, 
    0.015359942596352098`, 
    0.02448314785173377`}, {0.07381608109139871`, 
    0.05242856376555134`, 0.3026715722222353`, -0.0819524814567561`, 
    0.07919970253622352`, -0.033627335100296744`}, 
{0.01761648518373638`, -0.009970312195220557`, -0.08195248145675606`, 
    0.5519153257149109`, -0.1068423071796184`, 
    0.35028088018797104`}, {0.05925859721052879`, 
    0.015359942596352072`, 
    0.07919970253622348`, -0.10684230717961839`, 0.29853333010040306`,
     0.002543954098042138`}, {0.0458989604009327`, 
    0.02448314785173373`, -0.033627335100296814`, 
    0.35028088018797104`, 0.002543954098042013`, 0.6962012633235793`}};
ii = IdentityMatrix[d];
RandomVariate@MultinormalDistribution[sigma] (* fails *)
RandomVariate@MultinormalDistribution[sigma + ii] (* works *)

The example above was generated by code below, it’s a random covariance matrix with geometrically decaying eigenvalues, so not even that close to singular —

randomRotation[n_] := Module[{M, z, q, r, d, ph},
   z = RandomVariate[NormalDistribution[0, 1], {n, n}];
   {q, r} = QRDecomposition[z];
   d = Diagonal[r];
   ph = d/Abs[d];
   M = q*ph; (* Note, determinant may be -1 giving "pseudo-rotation", 
   switch 2 rows of the matrix to guarantee true rotation *)
   indices = If[Det[M] > 0, Range[n], {2, 1}~Join~Range[3, n]];
   M[[indices]]
   ];

d = 20;
rot = randomRotation[d];
ii = IdentityMatrix[d];

sigma = rot.DiagonalMatrix[Table[1/i, {i, 1, d}]].Inverse[rot];
RandomVariate@MultinormalDistribution[sigma] (* fails *)
RandomVariate@MultinormalDistribution[sigma+ii] (* works *)
```

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP