# Exsample of use of "nordpred"-functions. # Nordpred: R (www.r-project.org) & S-PLUS (www.insightful.com) functions # for prediction of cancer incidence (as used in the Nordpred project). # Written by: Bjørn Møller and Harald Fekjaer <harald.fekjar@kreftregisteret.no>, 2000-2002 # See also: http://www.kreftregisteret.no/software/nordpred # Reading package: source("Nordpred.S") # Reading data (Colon cancer for Norwegian males) indata <- read.table("data//colon-men-Norway.txt",header =T,sep=",",row.names=1) inpop1 <- read.table("data//men-Norway.txt",header =T,sep=",",row.names=1) inpop2 <- read.table("data//men-Norway-pred.txt",header =T,sep=",",row.names=1) # Include possible population predictions inpop <- cbind(inpop1,inpop2) # Run predictions: est <- nordpred.estimate(cases=indata,pyr=inpop,noperiod=4,startestage=5) res <- nordpred.prediction(est,startuseage=6,cuttrend=c(0,.25,.5,.75,.75),recent=T) # This can also be done in one command: res <- nordpred(cases=indata,pyr=inpop,startestage=5,startuseage=6,noperiods=4,cuttrend=c(0,.25,.5,.75,.75)) # The "nordpred"-function can also choose number periods to base predions on: # This is done by listing candidate number of periods in "noperiods". # If the goodness of fit test is rejected based on the widest base, # the first period is exclude etc. res <- nordpred(indata,inpop,startestage=5,startuseage=6,noperiods=4:6,cuttrend=c(0,.25,.5,.75,.75)) # Or with poisson link function (instead of the powerlink as used in the nordpred predictions): est2 <- nordpred.estimate(indata,inpop,4,5,linkfunc="poisson") res2 <- nordpred.prediction(est2,startuseage=6,cuttrend=c(0,.25,.5,.75,.75),recent=T) # Get results: print.nordpred(res) nordpred.getpred(res) summary(res,printpred=F) # Get results with stanardiziotion: wstand <- c(0.12, 0.1, 0.09, 0.09, 0.08, 0.08, 0.06, 0.06, 0.06, 0.06,0.05, 0.04, 0.04, 0.03, 0.02, 0.01, 0.005, 0.005) round(nordpred.getpred(res,incidence=T,standpop=NULL),2) round(nordpred.getpred(res,incidence=T,standpop=wstand),2) # Plot results: plot(res,standpop=wstand) # Plot results with power5 and poisson links: plot(res2,standpop=wstand) plot(res,new=F,lty=c(1,2),standpop=wstand) # Different cut trend scenarios, using average drift (recent=F): plot(nordpred.prediction(est,startuseage=6,cuttrend=c(0,0,0,0,0),recent=F),standpop=wstand,new=T) plot(nordpred.prediction(est,startuseage=6,cuttrend=c(1,1,1,1,1),recent=F),standpop=wstand,new=F,lty=c(1,2)) plot(nordpred.prediction(est,startuseage=6,cuttrend=c(0,.25,.5,.75,.75),recent=F),standpop=wstand,new=F,lty=c(1,4))