Accepted Solutions
Actually since (P1*P2*Q3+Q1*P2*Q3) can be reduced to (P1+Q1)*P2*Q3 and P1+Q1 is just 1 you can define C1 = 1 and then you don't need separate definition for C3. This eliminates the need for second lag of the p values.
You also don't really need to second P value since it is just the compliment of the first.
data want ;
set have ;
by id day;
length visit c lagp lagc lag2c 8;
if first.id then call missing(of visit lag:);
visit+1;
select (visit);
when (1) c=1;
when (2) c=p + lagp - p*lagp ;
otherwise c=lag2c*lagp*(1-p) + lagc*p ;
end;
output;
lagp=p;
lag2c=lagc;
lagc=c;
retain lag: ;
run;
You can simpify even more if you just initialize the lagged values to 1.
data want ;
set have ;
by id day;
if first.id then do;
visit=0; lagp=1; lagc=1; lag2c=1;
end;
visit+1;
c=lag2c*lagp*(1-p) + lagc*p ;
output;
lagp=p;
lag2c=lagc;
lagc=c;
retain visit lagp lagc lag2c ;
run;
You posted way more data than we need to test. Too many observations and too many variables.
It looks like you have just ID, VISIT, DAY , PROB and your EXPECTED value for C. The other probability is just 1-PROB and the other variables are just lags of the two prob variables which are easy to calculate.
data have ;
infile cards dsd dlm='|' truncover ;
input id visit day prob expect ;
cards;
1|1|0|0.999999838|
1|2|2|0.999991138|1
1|3|7|0.987787081|0.999999892
1|4|16|0.237545783|0.990688183
1|5|29|0.079696412|0.297568507
1|6|57|0.001868728|0.079362825
1|7|91|8.13E-6|0.000556716
1|8|120|1.15E-7|6.46E-7
1|9|150|6.64E-10|6.43E-11
1|10|180|1.76E-12|4.29E-16
1|11|210|2.11E-15|1.13E-22
1|12|240|0|9.05E-31
;
proc print;
format _numeric_ best32.;
run;
You said the problem was recursive. So to do that we need to know two things.
1) How to calculate the initial value. So it looks like when visit=1 it is unknown and when visit=2 it is 1.
2) How to calculate Nth value given the previous values.
From your post it looks like formula is based on the values from the previous 2 values some how.
I think I have translated this into code now. So here is sample input data. We really just need the ID, DAY (to insure the order is right) and P value. I added in your expected C series.
data have ;
infile cards dsd dlm='|' truncover ;
input id day p expect ;
cards;
1|0|0.999999838|
1|2|0.999991138|1
1|7|0.987787081|0.999999892
1|16|0.237545783|0.990688183
1|29|0.079696412|0.297568507
1|57|0.001868728|0.079362825
1|91|8.13E-6|0.000556716
1|120|1.15E-7|6.46E-7
1|150|6.64E-10|6.43E-11
1|180|1.76E-12|4.29E-16
1|210|2.11E-15|1.13E-22
1|240|0|9.05E-31
;Now set that data and you can calculate the variables needed for the formulas from those values. Instead of trying to get LAG() function to work I just used RETAIN. I moved the OUTPUT to before I generate the lagged values so you can see on the output the values used to calculate C.
data want ;
set have ;
by id day;
length visit c lagp lag2p q lagq lag2q lagc lag2c 8;
if first.id then call missing(of visit lag:);
visit+1;
q=1-p;
select (visit);
when (1) c=.;
when (2) c=lagp*p + lagq*p + lagp*q ;
when (3) c=lag2p*lagp*p + lag2q*lagp*p + lag2p*lagp*q + lag2q*lagp*q + lag2p*lagq*p;
otherwise c=lag2c*lagp*q + lagc*p ;
end;
output;
lag2p=lagp;
lagp=p;
lag2q=lagq;
lagq=q;
lag2c=lagc;
lagc=c;
retain lag: ;
run;
proc print; run;
Actually since (P1*P2*Q3+Q1*P2*Q3) can be reduced to (P1+Q1)*P2*Q3 and P1+Q1 is just 1 you can define C1 = 1 and then you don't need separate definition for C3. This eliminates the need for second lag of the p values.
You also don't really need to second P value since it is just the compliment of the first.
data want ;
set have ;
by id day;
length visit c lagp lagc lag2c 8;
if first.id then call missing(of visit lag:);
visit+1;
select (visit);
when (1) c=1;
when (2) c=p + lagp - p*lagp ;
otherwise c=lag2c*lagp*(1-p) + lagc*p ;
end;
output;
lagp=p;
lag2c=lagc;
lagc=c;
retain lag: ;
run;You can simpify even more if you just initialize the lagged values to 1.
data want ;
set have ;
by id day;
if first.id then do;
visit=0; lagp=1; lagc=1; lag2c=1;
end;
visit+1;
c=lag2c*lagp*(1-p) + lagc*p ;
output;
lagp=p;
lag2c=lagc;
lagc=c;
retain visit lagp lagc lag2c ;
run;
