# UseCase4 : Multiple administration routes - Oral and IV dosing. Defined by ODEs.
# Warfarin population pharmacokinetic model for multiple dosing and administration routes
# Dosing regimen: intravenous infusion followed by single oral administration
# Structural model : 1 compartment model (CL,V, FORAL and TLAG)
#    Structural model defined by ODEs
#    Single dose only.
# Covariate model : Weight on V and CL 
# Variability model: 
#   - IIV on CL, V, ka, FORAL (logit) and TLAG
#   - Correlation between CL and V random variables
#   - Combined residual error model


warfarin_PK_IVPO_dat = dataObj {
	DECLARED_VARIABLES{
		GUT :: dosingTarget
		CENTRAL :: dosingTarget
		Y ::observation
	}

  DATA_INPUT_VARIABLES{
      ID: 	{ use is id }
      TIME: { use is idv }
      WT: 	{ use is covariate }
      AMT: 	{ use is amt, define={1 in CMT as GUT, 2 in CMT as CENTRAL} }
      RATE: { use is rate }
      CMT: 	{ use is cmt }
      DV: 	{ use is dv, variable=Y}
      logtWT: { use is covariate }
      MDV: { use is mdv }
   }# end DATA_INPUT_VARIABLES

	DATA_DERIVED_VARIABLES{
		# This is an idea of how to solve the dose time problem.  Like 'use is amt' we assume that the DT variable is
		# only assigned when AMT > 0.  The typing ensured that the attributes reference a column with the correct 'use'.
		DT : { use is doseTime, idvColumn=TIME, amtColumn=AMT }
	}

   SOURCE{
      srcfile : { file="warfarin_infusion_oral.csv",
      			  inputFormat is nonmemFormat }
   }# end SOURCE
} # end data object

warfarin_PK_IVPO_prior = priorObj {

	PRIOR_PARAMETERS{
		MU_CL=ln(0.1)
		MU_V=ln(8)
		MU_KA=ln(0.362)
		MU_TLAG=ln(1)
		MU_PPV_CL=ln(0.1)
		MU_PPV_V=ln(0.1)
		MU_PPV_KA=ln(0.1)
		MU_FORAL=ln(0.7)
		MU_PPV_FORAL=ln(0.1)
		MU_CORR_CL_V=ln(0.01)
		MU_RUV_PROP=ln(0.1)
		MU_RUV_ADD=ln(0.1)

		VAR_CL=ln(0.1^2)^2
		VAR_V=ln(8^2)^2
		VAR_KA=ln(0.362^2)^2
		VAR_TLAG=ln(1^2)^2
		VAR_PPV_CL=ln(0.1^2)^2
		VAR_PPV_V=ln(0.1^2)^2
		VAR_PPV_KA=ln(0.1^2)^2
		VAR_FORAL=ln(0.7^2)^2
		VAR_PPV_FORAL=ln(0.1^2)^2
		VAR_CORR_CL_V=ln(0.01^2)^2
		VAR_RUV_PROP=ln(0.1^2)^2
		VAR_RUV_ADD=ln(0.1^2)^2
	}
	
	PRIOR_VARIABLE_DEFINITION{
		BETA_CL_WT=0.75
		BETA_V_WT=1
		PPV_TLAG=0.1
		
		POP_CL ~ LogNormal2(meanLog=MU_CL, varLog = VAR_CL)
		POP_V ~ LogNormal2(meanLog=MU_V, varLog = VAR_V)
		POP_KA ~ LogNormal2(meanLog=MU_KA, varLog = VAR_KA)
		POP_TLAG ~ LogNormal2(meanLog=MU_TLAG, varLog = VAR_TLAG)
		POP_FORAL ~ LogNormal2(meanLog=MU_FORAL, varLog = VAR_FORAL)
		PPV_CL ~ LogNormal2(meanLog=MU_PPV_CL, varLog = VAR_PPV_CL)
		PPV_V ~ LogNormal2(meanLog=MU_PPV_V, varLog = VAR_PPV_V)
		PPV_KA ~ LogNormal2(meanLog=MU_PPV_KA, varLog = VAR_PPV_KA)
		PPV_FORAL ~ LogNormal2(meanLog=MU_PPV_FORAL, varLog = VAR_PPV_FORAL)
		CORR_CL_V ~ LogNormal2(meanLog=MU_CORR_CL_V, varLog = VAR_CORR_CL_V)
		RUV_PROP ~ LogNormal2(meanLog=MU_RUV_PROP, varLog = VAR_RUV_PROP)
		RUV_ADD ~ LogNormal2(meanLog=MU_RUV_ADD, varLog = VAR_RUV_ADD)
	}
}


warfarin_PK_IVPO_par = parObj {
	# By default a parameter is to be estimated if fix is omitted
 	STRUCTURAL {
		POP_CL : { value = 0.1, lo = 0.001 }
		POP_V : { value = 8, lo = 0.001 } 
		POP_KA : { value = 0.362, lo = 0.001 }
		POP_TLAG : { value=1, lo=0.001 }
		BETA_CL_WT : { value = 0.75, fix = true }
		BETA_V_WT : { value = 1, fix = true } 
		POP_FORAL : { value=0.7, lo=0.001 }
		} # end STRUCTURAL
	VARIABILITY {
		PPV_CL : { value = 0.1 }
		PPV_V : { value = 0.1 }
		PPV_KA : { value = 0.1 }
		PPV_TLAG : { value = 0.1, fix=true }
		PPV_FORAL : { value = 0.1 }
		CORR_CL_V : { value = 0.01 }
		RUV_PROP : { value = 0.1, lo = 0 }
		RUV_ADD : { value = 0.1, lo = 0.0001 } 	
		} # end VARIABILITY
} # end of parameter object 

warfarin_PK_IVPO_mdl = mdlObj {
	IDV{ T }

	COVARIATES{
		logtWT
	}

	VARIABILITY_LEVELS{
		ID : { level=2, type is parameter }
		DV : { level=1, type is observation }
	} 

	STRUCTURAL_PARAMETERS { 
		POP_CL
		POP_V
		POP_KA
		POP_TLAG
		POP_FORAL
		BETA_CL_WT
		BETA_V_WT
	} # end STRUCTURAL_PARAMETERS
	
	VARIABILITY_PARAMETERS {
		PPV_CL
		PPV_V
		CORR_CL_V
		PPV_KA
		PPV_TLAG
		PPV_FORAL
		RUV_PROP
		RUV_ADD
	} # end VARIABILITY_PARAMETERS 
	
	RANDOM_VARIABLE_DEFINITION(level=ID) {
		ETA_CL ~ Normal(mean = 0, sd = PPV_CL)
		ETA_V ~ Normal(mean = 0, sd = PPV_V)
		:: {type is correlation, rv1=ETA_CL, rv2=ETA_V, value=CORR_CL_V}
		ETA_KA ~ Normal(mean = 0, sd = PPV_KA)
		ETA_TLAG ~ Normal(mean = 0, sd = PPV_TLAG) # We define correlation here
		ETA_FORAL ~ Normal(mean = 0, sd = PPV_FORAL)
	} # end RANDOM_VARIABLE_DEFINITION 

	
	INDIVIDUAL_VARIABLES { # This maps to the "Type 3" individual parameter definition in PharmML
	    CL : {type is linear, trans is ln, pop = POP_CL, fixEff = {coeff=BETA_CL_WT, cov=logtWT} , ranEff = ETA_CL }
	    V : {type is linear, trans is ln, pop = POP_V, fixEff = {coeff=BETA_V_WT, cov=logtWT} , ranEff = ETA_V }
	    KA : {type is linear, trans is ln, pop = POP_KA, ranEff = ETA_KA }
	    TLAG : {type is linear, trans is ln, pop = POP_TLAG, ranEff = ETA_TLAG }
	    FORAL : {type is linear, trans is logit, pop = POP_FORAL, ranEff = ETA_FORAL }
	} # end INDIVIDUAL_VARIABLES

	MODEL_PREDICTION {
		DT
		## Here there is discontinuity in RATEIN which is not well handled with differential equations.
		## Ideally we would handle this by defining input mechanisms via COMPARTMENTS.
		## See UseCase4_1.mdl, UseCase4_2.mdl and UseCase4_3.mdl 
		## Note also the model below does not handle multiple doses since then DT = 0
		##   and the RATEIN is reset =0.
		DEQ{
			RATEIN = if(T-DT >= TLAG) then GUT * KA
					 else 0
			GUT : { deriv =(- RATEIN), init = 0, x0 = 0 }
			CENTRAL : { deriv =(RATEIN * FORAL - CL * CENTRAL / V), init = 0, x0 = 0 } 
		}
	    CC = CENTRAL / V 
	} # end MODEL_PREDICTION
	
	RANDOM_VARIABLE_DEFINITION(level=DV){
	    EPS_Y ~ Normal(mean = 0, var = 1) # This maps the standard error model in PharmML. The error model is predefined.
	}

	OBSERVATION {
	    Y : {type is combinedError1, additive = RUV_ADD, proportional = RUV_PROP, eps = EPS_Y, prediction = CC } 
	} # end OBSERVATION
} # end of model object

warfarin_PK_ODE_task = taskObj {
	ESTIMATE{ 
		set algo is mcmc
        TARGET_SETTINGS(target="BUGS"){
            set 
            nchains = 1,
            burnin = 1,
            niter = 100,
        #    parameters = "LAMBDA0_POP,LAMBDA1_POP,K1_POP,K2_POP,W0_POP,LAMBDA0,K2,CV",
         #    thinupdater=1,
          #  thinsamples=5, # uno dei due
            odesolver="LSODA", #oppure LSODA
        #    dic="true",
        	winbugsgui="true" # oppure no
        }
}
} # end of task object

warfarin_PK_ODE_mog = mogObj {
    OBJECTS{
		warfarin_PK_IVPO_dat : { type is dataObj }
		warfarin_PK_IVPO_mdl : { type is mdlObj }
		warfarin_PK_IVPO_prior : { type is priorObj }
		warfarin_PK_ODE_task : { type is taskObj }
	}
}
