data(phoneframes)
A_in_frames<-cbind(1, DatA$Domain=="ab")
B_in_frames<-cbind(DatB$Domain=="ba",1)
Bdes_pps<-svydesign(id=~1, fpc=~ProbB, data=DatB,pps=ppsmat(PiklB))
Ades_pps <-svydesign(id=~1, fpc=~ProbA,data=DatA,pps=ppsmat(PiklA))
## optimal constant (Hartley) weighting
mf_pps<-multiframe(list(Ades_pps,Bdes_pps),list(A_in_frames,B_in_frames),theta=0.74)
svytotal(~Lei,mf_pps)
Awts<-cbind(1/DatA$ProbA, ifelse(DatA$ProbB==0,0,1/DatA$ProbB))
Bwts<-cbind(ifelse(DatB$ProbA==0,0,1/DatB$ProbA),1/DatB$ProbB )
## dividing by the expected number of selections (BKA or HH estimator)
mf_pps2<-multiframe(list(Ades_pps,Bdes_pps),list(Awts,Bwts),estimator="expected")
svymean(~Lei,mf_pps2)
## Metcalf and Scott approximation
DatB$Stratum<-10
DatB$Frame<-2
DatA$Frame<-1
Dat_both<-rbind(DatA,DatB)
frame_weights<-c(0.742,1-0.742)
Dat_both$fweights<-with(Dat_both, ifelse(Frame==1,
ifelse(Domain=="ab", frame_weights[1]*1/ProbA,1/ProbA),
ifelse(Domain=="ba", frame_weights[2]*1/ProbB, 1/ProbB)))
MSdesign<-svydesign(id=~1, strata=~Stratum, weights=~fweights,data=Dat_both)
svymean(~Lei,MSdesign)
Run the code above in your browser using DataLab