Search code examples
probabilitymaple

Where are these negatives coming from in Maple execution?


I am interested in simulating the phenomenon of "regression to the mean". Say a 0-1 vector V of length N is "gifted" if the number of 1s in V is greater than N/2 + 5*sqrt(N).

I want Maple to evaluate a string of M 0-1 lists, each of length N, to determine whether they are gifted. Then, given that list V[i] is gifted, I want to evaluate the probability that list V[i+1] is gifted.

So far my code is failing in a strange way. So far all the code is supposed to do is create the list of sums (called 'total') and the list 'g' which carries a 0 if total[i] <= N/2 + 5sqrt(N), and a 1 otherwise.

Here is the code:

RS:=proc(N) local ra,i:

ra:=rand(0..1):

[seq(ra(),i=1..N)]:

end:

Gift:=proc(N,M) local total, i, g :

total:=[seq(add(RS(N)),i=1..M)]:
g:=[seq(0,i=1..M)]:

for i from 1 to M do

if total[i] > (N/2 + 5*(N^(1/2))) then
g[i]:=1

fi:
od:

print(total, g)

end:

The trouble is, Maple responds, when I try Gift(100,20), "Error, (in Gift) cannot determine if this expression is true or false: 5*100^(1/2) < -2" or, when I try Gift(10000,20), "Error, (in Gift) cannot determine if this expression is true or false: 5*10000^(1/2) < -103."

Where are these negative numbers coming from? And why can't Maple tell whether 5(10000)^{1/2} < -103 or not?


Solution

  • The negative quantities are simply the part of the inequality that results when the portion with the radical is moved to one side and the purely rational portion is moved to the other.

    Use an appropriate mechanism for the resolution of the conditional test. For example,

    if is( total[i] > (N/2 + 5*N^(1/2)) ) then
      ...etc
    

    or, say,

    temp := evalf(N/2 + 5*N^(1/2));
    for i from 1 to M do
      if total[i] > temp then
        ...etc