41 enum stage_t {DiseaseFree,Precursor,PreClinical,Clinical,
Death};
44 enum event_t {toPrecursor, toPreClinical, toClinical, toDeath, Count};
47 string stage_names[5] = {
"DiseaseFree",
"Precursor",
"PreClinical",
"Clinical",
"Death"};
58 double Lam1,sigm1,p2,lam2,mu3,tau3;
63 static std::map<std::string, std::vector<double> >
report;
65 static void resetPopulation ();
69 CalibPerson(
double *par,
int i=0) {
81 virtual void handleMessage(
const cMessage* msg);
82 virtual Time age() {
return now(); }
85 void CalibPerson::resetPopulation() {
95 void CalibPerson::init() {
96 if (R::runif(0,1)<p2){
103 double lam1 = exp(R::rnorm(Lam1,sigm1));
104 scheduleAt(R::rexp(lam1), toPrecursor);
105 double x = R::runif(0,1);
106 scheduleAt((65 - 15*log(-log(x))), toDeath);
107 for(
int i=10;i<110;i=i+10){
115 void CalibPerson::handleMessage(
const cMessage* msg) {
117 double ctime[] = {20,40,60,80};
121 if (msg->
kind == toDeath) {
123 clinTime=std::min(clinTime,
now());
125 for(
unsigned int i=0; i<4 ; i++){
126 if(i <
report[
"TimeAtRisk"].size()){
127 report[
"TimeAtRisk"][i] += std::min(ctime[i],clinTime);
130 report[
"TimeAtRisk"].push_back(std::min(ctime[i],clinTime));
133 if(clinTime < ctime[i]){
142 else if (msg->
kind == toPrecursor) {
146 scheduleAt(dwellTime, toPreClinical);
150 else if (msg->
kind == toPreClinical) {
152 simtime_t dwellTime =
now()+ exp(R::rnorm(mu3,tau3*mu3));
153 scheduleAt(dwellTime, toClinical);
156 else if (msg->
kind == toClinical) {
159 string stagestr = stage_names[stage];
162 else if (msg->
kind == Count){
163 cind = min(9,
int(
now()/10 - 1));
164 string stagestr = stage_names[stage];
167 report[stagestr].assign(10,0);
169 report[stagestr][cind]+=1;
177 Rcpp::List parmsl(parms);
178 int nin = Rcpp::as<int>(parmsl[
"n"]);
179 std::vector<double> par = Rcpp::as<std::vector<double> >(parmsl[
"runpar"]);
181 CalibPerson::resetPopulation();
188 for (
int i = 0; i < nin; i++) {
189 person = CalibPerson(&par[0],0);