https://github.com/cran/rstpm2
Revision 643c36ce2deda0fb6cb6da4ee9fcefa95c8bdfa1 authored by Mark Clements on 26 July 2015, 18:57:30 UTC, committed by cran-robot on 26 July 2015, 18:57:30 UTC
0 parent
Raw File
Tip revision: 643c36ce2deda0fb6cb6da4ee9fcefa95c8bdfa1 authored by Mark Clements on 26 July 2015, 18:57:30 UTC
version 1.2.2
Tip revision: 643c36c
test.do

* stata
*cd c:\Users\marcle\Documents\Home\
clear
webuse brcancer
*use brcancer
stset rectime, f(censrec==1)
cap program drop dopredictions
program define dopredictions
  preserve
  predict hr, hrnumerator(hormon 1) ci
  predict haz, hazard ci
  predict surv, surv ci
  predict sdiff, sdiff1(hormon 1) ci
  list hr* in 1/5
  list haz* surv* in 1/5
  list sdiff* in 1/5
  restore
end

* basic model
stpm2 hormon, df(3) scale(h)
dopredictions

* logit model
stpm2 hormon, df(3) scale(odds)
dopredictions

* normal model
stpm2 hormon, df(3) scale(normal)
dopredictions

* tvc
stpm2 hormon, df(3) tvc(hormon) dftvc(3) scale(h)
dopredictions

* delayed entry
preserve
  replace _t0 = rectime*0.5 if hormon==0
  stpm2 hormon, df(3) scale(h)
  dopredictions
restore

* delayed entry and tvc
preserve
  replace _t0 = rectime*0.5 if hormon==0
  stpm2 hormon, df(3) scale(h) tvc(hormon) dftvc(3)
  dopredictions
restore

* cure
stpm2 hormon, df(3) scale(h) cure
dopredictions


* relative survival
preserve  
  gen rate0=10^(-5+x1/100)
  stpm2 hormon, df(3) scale(h) bhazard(rate0)
  dopredictions
restore

* test speed
clear all
set mem 100m
webuse brcancer
stset rectime, f(censrec==1)
expand 500
timer clear
timer on 1
stpm2 hormon, df(3) scale(h)
timer off 1
timer list

** Cure
* http://www.pauldickman.com/survival/solutions/q37.do
clear
use "http://www.pauldickman.com/survival/colon.dta"
stset surv_mm, failure(status=1 2) scale(12) exit(time 120.5)
gen _age = min(int(age + _t),99)
gen _year = int(yydx + _t)
sort _year sex _age
merge m:1 _year sex _age using "http://www.pauldickman.com/survival/popmort.dta",  keep(match master)

** Cure model using stpm2 **
stpm2 year8594, df(6) bhazard(rate) scale(hazard) cure 
predict surv, surv ci
list surv surv_lci surv_uci in 1/6

stpm2 year8594, df(6) bhazard(rate) scale(hazard)  
stpm2 year8594, df(6) scale(hazard)  

preserve
clear
set obs 101
gen x=_n-1
gen y=(x-50)^2
replace y=0 if _n==20 | _n==30
rcsgen x, knots(0 25 50 75 95 100) gen(rcs) reverse
return list
quietly su y
di `r(sum)'
regress y rcs*
quietly predict yhat
scatter x yhat 
restore
back to top