[R] r programming help III

From: Mohammad Ehsanul Karim <wildscop_at_yahoo.com>
Date: Fri 24 Jun 2005 - 19:15:56 EST


Dear List,

Is there any way we can make the following formula any shorter by means of R programming:

# Transition Probabilities observed for 19 years p01<-0.4052863; p11<-0.7309942

# Now we try to find expected probabilities
# for number of wet days in a week
# as defined by Gabriel, Neumann (1962)
N<-7

S<-1; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p1N0
S<-2; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p2N0
S<-3; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p3N0
S<-4; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p4N0
S<-5; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p5N0
S<-6; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p6N0
S<-7; c0<-N+1/2-abs(2*S-N-1/2);c<-1:c0;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose((S-1),(b-1))*choose((N-S),a)*(((1-p11)/(1-p01))^a)*((p01/p11)^b))->p7N0
# we obtain probabilities for dry case
pN0<-c(p11^N,p1N0,p2N0,p3N0,p4N0,p5N0,p6N0,p7N0)

S<-0; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p0N1
S<-1; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p1N1
S<-2; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p2N1
S<-3; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p3N1
S<-4; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p4N1
S<-5; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p5N1
S<-6; c1<-N+1/2-abs(2*S-N+1/2);c<-1:c1;
a<-ceiling((c-1)/2); b<-ceiling(c/2);
(p11^S)*((1-p01)^(N-S))*sum(choose(S,a)*choose((N-S-1),(b-1))*(((1-p11)/(1-p01))^b)*((p01/p11)^a))->p6N1
# we obtain probabilities for wet case
pN1<-c(p0N1,p1N1,p2N1,p3N1,p4N1,p5N1,p6N1, p11^N) # Use of transition probabilities to get expected probabilities
d<-p11-p01
P<-p01/(1-d)
pN<-(1-P)*pN0+P*pN1
# Expected Probabilities for number of wet days in a week is pN
pN # gives the following output:

[1] 0.05164856 0.05126159 0.10424380 0.16317247 0.20470403 0.20621178 0.16104972
[8] 0.09170593

Any hint, help, support, references will be highly appreciated.
Thank you for your time.


Mohammad Ehsanul Karim

Web: http://snipurl.com/ehsan
ISRT, University of Dhaka, BD



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Fri Jun 24 19:22:39 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:33:02 EST