Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Working on some large systems using DiscreteMarkovProcess, I changed the transition matrix to machine precision vs using exact values, which sped things up handily. The only problem was for edge cases, accuracy suffered. I expected that. However, when probing the limits using simple cases, I was surprised to find this happens when I'd not expected it.

A simplified test example:

Table[
 intleg = tleg/2;

 sa = SparseArray[{{i_, i_} :> 
     1 - (intleg - i + 1)/tleg, {i_, j_} /; 
      j == i + 1 :> (1 - i + intleg)/tleg}, {intleg + 1, intleg + 1}];


 saN = N[sa, 50];

 mp = DiscreteMarkovProcess[1, sa];
 mpN = DiscreteMarkovProcess[1, saN];

 {tleg, N@Mean[FirstPassageTimeDistribution[mp, intleg + 1]],
  N@Mean[FirstPassageTimeDistribution[mpN, intleg + 1]]},

 {tleg, {10, 20, 30, 40}}]

(* {{10,22.8333,22.8333},{20,58.5794,58.5793},{30,99.5469,-1.24192},{40,143.91,0.}} *)

This simulates a simple "coupon-collector" situation with a population of tleg distinct items, where only half are "interesting" to me. So a transition matrix is built with appropriate probabilities with an absorbing state when all "interesting" are obtained. That same matrix is converted to inexact via N.

Comparing the results, we see the last two are not even in the zip-code of the correct result. Upping the accuracy/precision of the conversion via N seems to have no effect.

Any ideas what's going on here (I suspect some kind of catastrophic cancellation because of the characteristics of the probabilities)?

share|improve this question
    
I really don't know. It looks like it could be a bug rather than a numerical problem, as the Precision of the results doesn't decrease very much with respect to the input. But possibly this is just a case that is handled very badly by significance arithmetic. I suspect we need @AndyRoss to tell us what's really going on here? –  Oleksandr R. Jun 18 at 0:16
1  
This one wasn't my baby so I can't really say for sure what is going on. I suspect that a completely different algorithm is being used for inexact numbers and it is failing. You should file a bug. –  Andy Ross Jun 18 at 1:54
    
@AndyRoss: I shall do so, however I'll wait a bit to see and responses and if perhaps there's something I missed... Thanks for reply. –  rasher Jun 18 at 4:56
    
@OleksandrR.: Thanks for reply - yes, it is a strange one, I'm still poking at it, but see AndyRoss' reply... –  rasher Jun 18 at 4:56

1 Answer 1

up vote 4 down vote accepted

Thanks for the report, @rasher. This was a bug and it is already fixed in the development version.

In[1]:= Table[intleg = tleg/2;
 sa = SparseArray[{{i_, i_} :> 
     1 - (intleg - i + 1)/tleg, {i_, j_} /; 
      j == i + 1 :> (1 - i + intleg)/tleg}, {intleg + 1, intleg + 1}];
 saN = N[sa, 50];
 mp = DiscreteMarkovProcess[1, sa];
 mpN = DiscreteMarkovProcess[1, saN];
 {tleg, N@Mean[FirstPassageTimeDistribution[mp, intleg + 1]], 
  N@Mean[FirstPassageTimeDistribution[mpN, intleg + 1]]}, {tleg, {10, 20, 30, 40}}]

Out[1]= {{10, 22.8333, 22.8333}, {20, 58.5794, 58.5794}, {30, 99.5469, 99.5469}, {40, 143.91, 143.91}}

Unfortunately I don't have a workaround at the moment. I apologize for the inconvenience.

share|improve this answer
    
Thanks for confirmation (Are you WRI staff?). Glad to hear it's fixed in upcoming release. –  rasher Jun 19 at 21:07
    
Yes, I worked on the Markov process functionality. –  Bhuvanesh Jun 19 at 21:33
    
Cool! Thanks again for response, I'll go ahead and add bug tag. –  rasher Jun 19 at 21:34

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.