[R] r programming help III

Mohammad Ehsanul Karim wildscop at yahoo.com
Fri Jun 24 11:15:56 CEST 2005


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




More information about the R-help mailing list