Index: tol/stdlib/math/stat/models/bayesian/bysMcmc/bsr/gibbs/_basicMaster.tol
===================================================================
--- tol/stdlib/math/stat/models/bayesian/bysMcmc/bsr/gibbs/_basicMaster.tol	(revision d6a7fad92a516aed1f29b5d544ea3055eb287d4a)
+++ tol/stdlib/math/stat/models/bayesian/bysMcmc/bsr/gibbs/_basicMaster.tol	(revision b925b0b2a6013e16d99e2614d8c28d1905d46c7e)
@@ -459,5 +459,6 @@
       Matrix filter = IfMat(filter.mis_,RProd(filter_,?),filter_);
       Matrix filter.mis  = IfMat(filter.mis_,filter.mis_,RProd(filter.mis_,?));
-      Matrix noise = VMat2Mat(SubRow(fullNoise,segment->EquIdx));
+      VMatrix noise.v = SubRow(fullNoise,segment->EquIdx);
+      Matrix noise = VMat2Mat(noise.v);
     //WriteLn("TRACE"+_MID+" eval 7["<<k+"].2.2 segment->SigmaIndex="<<segment->SigmaIndex);
       Real L.det = 
@@ -465,7 +466,7 @@
         //The determinant of a triangular matrix is the product of its diagonal
         VMatrix L.diag = SubRow(Tra(SubDiag(segment->L,0)),segment->EquIdx);
-        Exp(VMatSum(Log(L.diag)))^(1/Card(segment->EquIdx))
+        Exp(VMatAvr(Log(L.diag)))
       };
-      Real sigma  = L.det*If(segment->SigmaIndex<=0, 1,
+      Real sigma  = If(segment->SigmaIndex<=0, 1,
        Sqrt(VMatDat(_.sig.blk::_.store,segment->SigmaIndex,1)));
     //WriteLn("TRACE"+_MID+" eval 7["<<k+"].2.3");
@@ -476,5 +477,6 @@
 
       Matrix difNoise = noise;
-      Matrix residuals = noise;
+      VMatrix Li = SubRow(SubCol(segment->Li,segment->EquIdx),segment->EquIdx);
+      Matrix residuals = VMat2Mat(Li * noise.v);
       Matrix forecast = output.full - noise;
       Real difDeg = 0;