[R] DeSolver giving "NA" as output, but running fully.

John Kane jrkrideau at inbox.com
Mon Apr 27 23:46:46 CEST 2015


Data?

Use deput()  (see ?dput) to provide some sample data.  Also you might find this useful http://adv-r.had.co.nz/Reproducibility.html 

John Kane
Kingston ON Canada


> -----Original Message-----
> From: walke554 at umn.edu
> Sent: Mon, 27 Apr 2015 13:34:54 -0700 (PDT)
> To: r-help at r-project.org
> Subject: [R] DeSolver giving "NA" as output, but running fully.
> 
> Hello All,
> 
> I am currently looking on a transmission model for STD transmission
> within a
> population.  I am able to run my full code and the ODE function, but when
> I
> look at my output, all I get is "NA" for each time step beyond the first.
> There doesn't seem to be any syntax error, and I do get my entire program
> to
> run.  Here is my code:
> 
> setwd("C:/Users/L/Documents/MastersThesis")
> 
> require(deSolve);
> 
> ########
> #Model 1
> ########
> 
> #The function
> HPVInfection<-function(t,y,p){
> 	XFL = y[1]; #number of susceptible unvaccinated females low risk
> 	XFM = y[2]; #number of susceptible unvaccinated females medium risk
> 	XFH = y[3]; #number of susceptible unvaccinated females high risk
> 	XML = y[4]; #number of susceptible unvaccinated males low risk
> 	XMM = y[5]; #number of susceptible unvaccinated males medium risk
> 	XMH = y[6]; #number of susceptible unvaccinated males high risk
> 	Y1FL = y[7]; #number of infected unvaccinated females low risk infected
> with vaccine strain
> 	Y1FM = y[8]; #number of infected unvaccinated females medium risk low
> risk
> infected with vaccine strain
> 	Y1FH = y[9]; #number of infected unvaccinated females high risk low risk
> infected with vaccine strain
> 	Y1ML = y[10]; #number of infected unvaccinated males low risk low risk
> infected with vaccine strain
> 	Y1MM = y[11]; #number of infected unvaccinated males medium risk low
> risk
> infected with vaccine strain
> 	Y1MH = y[12]; #number of infected unvaccinated males high risk low risk
> infected with vaccine strain
> 	Y2FL = y[13]; #number of infected unvaccinated females low risk infected
> with non-vaccine strain
> 	Y2FM = y[14]; #number of infected unvaccinated females medium risk low
> risk
> infected with non-vaccine strain
> 	Y2FH = y[15]; #number of infected unvaccinated females high risk low
> risk
> infected with non-vaccine strain
> 	Y2ML = y[16]; #number of infected unvaccinated males low risk low risk
> infected with non-vaccine strain
> 	Y2MM = y[17]; #number of infected unvaccinated males medium risk low
> risk
> infected with non-vaccine strain
> 	Y2MH = y[18]; #number of infected unvaccinated males high risk low risk
> infected with non-vaccine strain
> 	ZFL = y[19]; #number of immune females low risk
> 	ZFM = y[20]; #number of immune females medium risk
> 	ZFH = y[21]; #number of immune females high risk
> 	ZML = y[22]; #number of immune males low risk
> 	ZMM = y[23]; #number of immune males medium risk
> 	ZMH = y[24]; #number of immune males high risk
> 	VFL = y[25]; #number of susceptible vaccinated females low risk
> 	VFM = y[26]; #number of susceptible vaccinated females medium risk
> 	VFH = y[27]; #number of susceptible vaccinated females high risk
> 	VML = y[28]; #number of susceptible vaccinated males low risk
> 	VMM = y[29]; #number of susceptible vaccinated males medium risk
> 	VMH = y[30]; #number of susceptible vaccinated males high risk
> 	W1FL = y[31]; #number of infected vaccinated females low risk infected
> with
> vaccine strain
> 	W1FM = y[32]; #number of infected vaccinated females medium risk
> infected
> with vaccine strain
> 	W1FH = y[33]; #number of infected vaccinated females high risk infected
> with vaccine strain
> 	W1ML = y[34]; #number of infected vaccinated males low risk infected
> with
> vaccine strain
> 	W1MM = y[35]; #number of infected vaccinated males medium risk infected
> with vaccine strain
> 	W1MH = y[36]; #number of infected vaccinated males high risk infected
> with
> vaccine strain
> 	W2FL = y[37]; #number of infected vaccinated females low risk infected
> with
> non-vaccine strain
> 	W2FM = y[39]; #number of infected vaccinated females medium risk
> infected
> with non-vaccine strain
> 	W2FH = y[40]; #number of infected vaccinated females high risk infected
> with non-vaccine strain
> 	W2ML = y[41]; #number of infected vaccinated males low risk infected
> with
> non-vaccine strain
> 	W2MM = y[42]; #number of infected vaccinated males medium risk infected
> with non-vaccine strain
> 	W2MH = y[43]; #number of infected vaccinated males high risk infected
> with
> non-vaccine strain
> 	with(as.list(p), {
> 		dXFL.dt = (0.5 * mew * omega[1,1] * (1-phi) * total) - ((partner[1,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (tau[1,1] * W1ML) + (tau[1,2] * W2ML)) /
> population[1,2]) * rho[1,1]) + (((Y1MM + Y2MM + (tau[1,1]*W1MM) +
> (tau[1,2]
> * W2MM))/population[2,2]) * rho[1,2]) + (((Y1MH + Y2MH + (tau[1,2] *
> W1MH) +
> (tau[1,2] * W2MH)) / population[3,2]) * rho[1,3])) + mew) * XFL) + (sigma
> *
> VFL);
> 		dXFM.dt = (0.5 * mew * omega[2,1] * (1-phi) * total) - ((partner[2,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (tau[1,1] * W1ML) + (tau[1,2] * W2ML)) /
> population[1,2]) * rho[2,1]) + (((Y1MM + Y2MM + (tau[1,1]*W1MM) +
> (tau[1,2]
> * W2MM))/population[2,2]) * rho[2,2]) + (((Y1MH + Y2MH + (tau[1,1] *
> W1MH) +
> (tau[1,2] * W2MH)) / population[3,2]) * rho[2,3])) + mew) * XFM) + (sigma
> *
> VFM);
> 		dXFH.dt = (0.5 * mew * omega[3,1] * (1-phi) * total) - ((partner[3,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (tau[1,1] * W1ML) + (tau[1,2] * W2ML)) /
> population[1,2]) * rho[3,1]) + (((Y1MM + Y2MM + (tau[1,1]*W1MM) +
> (tau[1,2]
> * W2MM))/population[2,2]) * rho[3,2]) + (((Y1MH + Y2MH + (tau[1,1] *
> W1MH) +
> (tau[1,2] * W2MH)) / population[3,2]) * rho[3,3])) + mew) * XFH) + (sigma
> *
> VFH);
> 		dXML.dt = (0.5 * mew * omega[1,1] * (1-phi) * total) - ((partner[1,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (tau[2,1] * W1FL) + (tau[2,2] * W2FL)) /
> population[1,1]) * rho[1,1]) + (((Y1FM + Y2FM + (tau[2,1]*W1FM) +
> (tau[2,2]
> * W2FM))/population[2,1]) * rho[1,2]) + (((Y1FH + Y2FH + (tau[2,1] *
> W1FH) +
> (tau[2,2] * W2FH)) / population[3,1]) * rho[1,3])) + mew) * XML) + (sigma
> *
> VML);
> 		dXMM.dt = (0.5 * mew * omega[2,1] * (1-phi) * total) - ((partner[2,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (tau[2,1] * W1FL) + (tau[2,2] * W2FL)) /
> population[1,1]) * rho[2,1]) + (((Y1FM + Y2FM + (tau[2,1]*W1FM) +
> (tau[2,2]
> * W2FM))/population[2,1]) * rho[2,2]) + (((Y1FH + Y2FH + (tau[2,1] *
> W1FH) +
> (tau[2,2] * W2FH)) / population[3,1]) * rho[2,3])) + mew) * XMM) + (sigma
> *
> VMM);
> 		dXMH.dt = (0.5 * mew * omega[3,1] * (1-phi) * total) - ((partner[3,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (tau[2,1] * W1FL) + (tau[2,2] * W2FL)) /
> population[1,1]) * rho[3,1]) + (((Y1FM + Y2FM + (tau[2,1]*W1FM) +
> (tau[2,2]
> * W2FM))/population[2,1]) * rho[3,2]) + (((Y1FH + Y2FH + (tau[2,1] *
> W1FH) +
> (tau[2,2] * W2FH)) / population[3,1]) * rho[3,3])) + mew) * XMH) + (sigma
> *
> VMH);
> 		dY1FL.dt = (XFL * (partner[1,1] * beta[1,1] * ((rho[1,1] * ((Y1ML +
> (tau[1,1] * W1ML)) / population[1,2])) + (rho[1,2] * ((Y1MM + (tau[1,1] *
> W1MM)) / population[2,2])) + (rho[1,3] * ((Y1MH + (tau[1,1] * W1MH)) /
> population[3,2]))))) - ((mew + gamma[2,1]) * Y1FL);
> 		dY1FM.dt = (XFM * (partner[2,1] * beta[1,1] * ((rho[2,1] * ((Y1ML +
> (tau[1,1] * W1ML)) / population[1,2])) + (rho[2,2] * ((Y1MM + (tau[1,1] *
> W1MM)) / population[2,2])) + (rho[2,3] * ((Y1MH + (tau[1,1] * W1MH)) /
> population[3,2]))))) - ((mew + gamma[2,1]) * Y1FM);
> 		dY1FH.dt = (XFH * (partner[3,1] * beta[1,1] * ((rho[3,1] * ((Y1ML +
> (tau[1,1] * W1ML)) / population[1,2])) + (rho[2,2] * ((Y1MM + (tau[1,1] *
> W1MM)) / population[2,2])) + (rho[3,3] * ((Y1MH + (tau[1,1] * W1MH)) /
> population[3,2]))))) - ((mew + gamma[2,1]) * Y1FH);
> 		dY1ML.dt = (XML * (partner[1,1] * beta[2,1] * ((rho[1,1] * ((Y1FL +
> (tau[2,1] * W1FL)) / population[1,1])) + (rho[1,2] * ((Y1FM + (tau[2,1] *
> W1FM)) / population[2,1])) + (rho[1,3] * ((Y1FH + (tau[2,1] * W1FH)) /
> population[3,1]))))) - ((mew + gamma[1,1]) * Y1MM);
> 		dY1MM.dt = (XMM * (partner[2,1] * beta[2,1] * ((rho[2,1] * ((Y1FL +
> (tau[2,1] * W1FL)) / population[1,1])) + (rho[2,2] * ((Y1FM + (tau[2,1] *
> W1FM)) / population[2,1])) + (rho[2,3] * ((Y1FH + (tau[2,1] * W1FH)) /
> population[3,1]))))) - ((mew + gamma[1,1]) * Y1MM);
> 		dY1MH.dt = (XMH * (partner[3,1] * beta[2,1] * ((rho[3,1] * ((Y1FL +
> (tau[2,1] * W1FL)) / population[1,1])) + (rho[2,2] * ((Y1FM + (tau[2,1] *
> W1FM)) / population[2,1])) + (rho[3,3] * ((Y1FH + (tau[2,1] * W1MH)) /
> population[3,1]))))) - ((mew + gamma[1,1]) * Y1MH);
> 		dY2FL.dt = (XFL * (partner[1,1] * beta[1,1] * ((rho[1,1] * ((Y2ML +
> (tau[1,2] * W2ML)) / population[1,2])) + (rho[1,2] * ((Y1MM + (tau[1,2] *
> W2MM)) / population[2,2])) + (rho[1,3] * ((Y2MH + (tau[1,2] * W2MH)) /
> population[3,2]))))) - ((mew + gamma[2,2]) * Y1FL);
> 		dY2FM.dt = (XFM * (partner[2,1] * beta[1,1] * ((rho[2,1] * ((Y2ML +
> (tau[1,2] * W2ML)) / population[1,2])) + (rho[2,2] * ((Y1MM + (tau[1,2] *
> W2MM)) / population[2,2])) + (rho[2,3] * ((Y2MH + (tau[1,2] * W2MH)) /
> population[3,2]))))) - ((mew + gamma[2,2]) * Y1FM);
> 		dY2FH.dt = (XFH * (partner[3,1] * beta[1,1] * ((rho[3,1] * ((Y2ML +
> (tau[1,2] * W2ML)) / population[1,2])) + (rho[2,2] * ((Y1MM + (tau[1,2] *
> W2MM)) / population[2,2])) + (rho[3,3] * ((Y2MH + (tau[1,2] * W2MH)) /
> population[3,2]))))) - ((mew + gamma[2,2]) * Y1FH);
> 		dY2ML.dt = (XML * (partner[1,1] * beta[2,1] * ((rho[1,1] * ((Y2FL +
> (tau[2,2] * W2FL)) / population[1,1])) + (rho[1,2] * ((Y1FM + (tau[2,2] *
> W2FM)) / population[2,1])) + (rho[1,3] * ((Y2FH + (tau[2,2] * W2FH)) /
> population[3,1]))))) - ((mew + gamma[1,2]) * Y1MM);
> 		dY2MM.dt = (XMM * (partner[2,1] * beta[2,1] * ((rho[2,1] * ((Y2FL +
> (tau[2,2] * W2FL)) / population[1,1])) + (rho[2,2] * ((Y1FM + (tau[2,2] *
> W2FM)) / population[2,1])) + (rho[2,3] * ((Y2FH + (tau[2,2] * W2FH)) /
> population[3,1]))))) - ((mew + gamma[1,2]) * Y1MM);
> 		dY2MH.dt = (XMH * (partner[3,1] * beta[2,1] * ((rho[3,1] * ((Y2FL +
> (tau[2,2] * W2FL)) / population[1,1])) + (rho[2,2] * ((Y1FM + (tau[2,2] *
> W2FM)) / population[2,1])) + (rho[3,3] * ((Y2FH + (tau[2,2] * W2MH)) /
> population[3,1]))))) - ((mew + gamma[1,2]) * Y1MH);
> 		dZFL.dt = ((gamma[2,1] * Y1FL) + (gamma[2,2] * Y2FL) + (gamma[4,1] *
> W1FL)
> + (gamma[4,2]*W2FL)) - (mew * ZFL);
> 		dZFM.dt = ((gamma[2,1] * Y1FM) + (gamma[2,2] * Y2FM) + (gamma[4,1] *
> W1FM)
> + (gamma[4,2]*W2FM)) - (mew * ZFM);
> 		dZFH.dt = ((gamma[2,1] * Y1FH) + (gamma[2,2] * Y2FH) + (gamma[4,1] *
> W1FH)
> + (gamma[4,2]*W2FH)) - (mew * ZFH);
> 		dZML.dt = ((gamma[1,1] * Y1ML) + (gamma[1,2] * Y2ML) + (gamma[3,1] *
> W1FM)
> + (gamma[3,2]*W2ML)) - (mew * ZML);
> 		dZMM.dt = ((gamma[1,1] * Y1MM) + (gamma[1,2] * Y2MM) + (gamma[3,1] *
> W1MM)
> + (gamma[3,2]*W2MM)) - (mew * ZMM);
> 		dZMH.dt = ((gamma[1,1] * Y1MH) + (gamma[1,2] * Y2MH) + (gamma[3,1] *
> W1MH)
> + (gamma[3,2]*W2MH)) - (mew * ZMH);
> 		dVFL.dt = (0.5 * mew * omega[1,1] * phi * total) - ((partner[1,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (delta[1,1] * tau[1,1] * W1ML) +
> (delta[1,2] *
> tau[1,2] * W2ML)) / population[1,2]) * rho[1,1]) + (((Y1MM + Y2MM +
> (delta[1,1] * tau[1,1]*W1MM) + (delta[1,2] * tau[1,2] *
> W2MM))/population[2,2]) * rho[1,2]) + (((Y1MH + Y2MH + (delta[1,1] *
> tau[1,2] * W1MH) + (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]) *
> rho[1,3])) + mew) * VFL) - (sigma * VFL);
> 		dVFM.dt = (0.5 * mew * omega[2,1] * phi * total) - ((partner[2,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (delta[1,1] * tau[1,1] * W1ML) +
> (delta[1,2] *
> tau[1,2] * W2ML)) / population[1,2]) * rho[2,1]) + (((Y1MM + Y2MM +
> (delta[1,1] * tau[1,1]*W1MM) + (delta[1,2] * tau[1,2] *
> W2MM))/population[2,2]) * rho[2,2]) + (((Y1MH + Y2MH + (delta[1,1] *
> tau[1,1] * W1MH) + (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]) *
> rho[2,3])) + mew) * VFM) - (sigma * VFM);
> 		dVFH.dt = (0.5 * mew * omega[3,1] * phi * total) - ((partner[3,1] *
> beta[1,1] * ((((Y1ML + Y2ML + (delta[1,1] * tau[1,1] * W1ML) +
> (delta[1,2] *
> tau[1,2] * W2ML)) / population[1,2]) * rho[3,1]) + (((Y1MM + Y2MM +
> (delta[1,1] * tau[1,1]*W1MM) + (delta[1,2] * tau[1,2] *
> W2MM))/population[2,2]) * rho[3,2]) + (((Y1MH + Y2MH + (delta[1,1] *
> tau[1,1] * W1MH) + (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]) *
> rho[3,3])) + mew) * VFH) - (sigma * VFH);
> 		dVML.dt = (0.5 * mew * omega[1,1] * phi * total) - ((partner[1,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (delta[2,1] * tau[2,1] * W1FL) +
> (delta[2,2] *
> tau[2,2] * W2FL)) / population[1,1]) * rho[1,1]) + (((Y1FM + Y2FM +
> (delta[2,1] * tau[2,1]*W1FM) + (delta[2,2] * tau[2,2] *
> W2FM))/population[2,1]) * rho[1,2]) + (((Y1FH + Y2FH + (delta[2,1] *
> tau[2,1] * W1FH) + (delta[2,2] * tau[2,2] * W2FH)) / population[3,1]) *
> rho[1,3])) + mew) * VML) - (sigma * VML);
> 		dVMM.dt = (0.5 * mew * omega[2,1] * phi * total) - ((partner[2,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (delta[2,1] * tau[2,1] * W1FL) +
> (delta[2,2] *
> tau[2,2] * W2FL)) / population[1,1]) * rho[2,1]) + (((Y1FM + Y2FM +
> (delta[2,1] * tau[2,1]*W1FM) + (delta[2,2] * tau[2,2] *
> W2FM))/population[2,1]) * rho[2,2]) + (((Y1FH + Y2FH + (delta[2,1] *
> tau[2,1] * W1FH) + (delta[2,2] * tau[2,2] * W2FH)) / population[3,1]) *
> rho[2,3])) + mew) * VMM) - (sigma * VMM);
> 		dVMH.dt = (0.5 * mew * omega[3,1] * phi * total) - ((partner[3,1] *
> beta[2,1] * ((((Y1FL + Y2FL + (delta[2,1] * tau[2,1] * W1FL) +
> (delta[2,2] *
> tau[2,2] * W2FL)) / population[1,1]) * rho[3,1]) + (((Y1FM + Y2FM +
> (delta[2,1] * tau[2,1]*W1FM) + (delta[2,2] * tau[2,2] *
> W2FM))/population[2,1]) * rho[3,2]) + (((Y1FH + Y2FH + (delta[2,1] *
> tau[2,1] * W1FH) + (delta[2,2] * tau[2,2] * W2FH)) / population[3,1]) *
> rho[3,3])) + mew) * VMH) - (sigma * VMH);
> 		dW1FL.dt = (VFL * (partner[1,1] * beta[1,1] * ((rho[1,1] * ((Y1ML +
> (delta[1,1] * tau[1,1] * W1ML)) / population[1,2])) + (rho[1,2] * ((Y1MM
> +
> (delta[1,1] * tau[1,1] * W1MM)) / population[2,2])) + (rho[1,3] * ((Y1MH
> +
> (delta[1,1] * tau[1,1] * W1MH)) / population[3,2]))))) - ((mew +
> gamma[4,1])
> * W1FL);
> 		dW1FM.dt = (VFM * (partner[2,1] * beta[1,1] * ((rho[2,1] * ((Y1ML +
> (delta[1,1] * tau[1,1] * W1ML)) / population[1,2])) + (rho[2,2] * ((Y1MM
> +
> (delta[1,1] * tau[1,1] * W1MM)) / population[2,2])) + (rho[2,3] * ((Y1MH
> +
> (delta[1,1] * tau[1,1] * W1MH)) / population[3,2]))))) - ((mew +
> gamma[4,1])
> * W1FM);
> 		dW1FH.dt = (VFH * (partner[3,1] * beta[1,1] * ((rho[3,1] * ((Y1ML +
> (delta[1,1] * tau[1,1] * W1ML)) / population[1,2])) + (rho[2,2] * ((Y1MM
> +
> (delta[1,1] * tau[1,1] * W1MM)) / population[2,2])) + (rho[3,3] * ((Y1MH
> +
> (delta[1,1] * tau[1,1] * W1MH)) / population[3,2]))))) - ((mew +
> gamma[4,1])
> * W1FH);
> 		dW1ML.dt = (VML * (partner[1,1] * beta[2,1] * ((rho[1,1] * ((Y1FL +
> (delta[2,1] * tau[2,1] * W1FL)) / population[1,1])) + (rho[1,2] * ((Y1FM
> +
> (delta[2,1] * tau[2,1] * W1FM)) / population[2,1])) + (rho[1,3] * ((Y1FH
> +
> (delta[2,1] * tau[2,1] * W1FH)) / population[3,1]))))) - ((mew +
> gamma[3,1])
> * W1MM);
> 		dW1MM.dt = (VMM * (partner[2,1] * beta[2,1] * ((rho[2,1] * ((Y1FL +
> (delta[2,1] * tau[2,1] * W1FL)) / population[1,1])) + (rho[2,2] * ((Y1FM
> +
> (delta[2,1] * tau[2,1] * W1FM)) / population[2,1])) + (rho[2,3] * ((Y1FH
> +
> (delta[2,1] * tau[2,1] * W1FH)) / population[3,1]))))) - ((mew +
> gamma[3,1])
> * W1MM);
> 		dW1MH.dt = (VMH * (partner[3,1] * beta[2,1] * ((rho[3,1] * ((Y1FL +
> (delta[2,1] * tau[2,1] * W1FL)) / population[1,1])) + (rho[2,2] * ((Y1FM
> +
> (delta[2,1] * tau[2,1] * W1FM)) / population[2,1])) + (rho[3,3] * ((Y1FH
> +
> (delta[2,1] * tau[2,1] * W1MH)) / population[3,1]))))) - ((mew +
> gamma[3,1])
> * W1MH);
> 		dW2FL.dt = (VFL * (partner[1,1] * beta[1,1] * ((rho[1,1] * ((Y2ML +
> (delta[1,2] * tau[1,2] * W2ML)) / population[1,2])) + (rho[1,2] * ((Y1MM
> +
> (delta[1,2] * tau[1,2] * W2MM)) / population[2,2])) + (rho[1,3] * ((Y2MH
> +
> (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]))))) - ((mew +
> gamma[4,2])
> * W1FL);
> 		dW2FM.dt = (VFM * (partner[2,1] * beta[1,1] * ((rho[2,1] * ((Y2ML +
> (delta[1,2] * tau[1,2] * W2ML)) / population[1,2])) + (rho[2,2] * ((Y1MM
> +
> (delta[1,2] * tau[1,2] * W2MM)) / population[2,2])) + (rho[2,3] * ((Y2MH
> +
> (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]))))) - ((mew +
> gamma[4,2])
> * W1FM);
> 		dW2FH.dt = (VFH * (partner[3,1] * beta[1,1] * ((rho[3,1] * ((Y2ML +
> (delta[1,2] * tau[1,2] * W2ML)) / population[1,2])) + (rho[2,2] * ((Y1MM
> +
> (delta[1,2] * tau[1,2] * W2MM)) / population[2,2])) + (rho[3,3] * ((Y2MH
> +
> (delta[1,2] * tau[1,2] * W2MH)) / population[3,2]))))) - ((mew +
> gamma[4,2])
> * W1FH);
> 		dW2ML.dt = (VML * (partner[1,1] * beta[2,1] * ((rho[1,1] * ((Y2FL +
> (delta[2,2] * tau[2,2] * W2FL)) / population[1,1])) + (rho[1,2] * ((Y1FM
> +
> (delta[2,2] * tau[2,2] * W2FM)) / population[2,1])) + (rho[1,3] * ((Y2FH
> +
> (delta[2,2] * tau[2,2] * W2FH)) / population[3,1]))))) - ((mew +
> gamma[3,2])
> * W1MM);
> 		dW2MM.dt = (VMM * (partner[2,1] * beta[2,1] * ((rho[2,1] * ((Y2FL +
> (delta[2,2] * tau[2,2] * W2FL)) / population[1,1])) + (rho[2,2] * ((Y1FM
> +
> (delta[2,2] * tau[2,2] * W2FM)) / population[2,1])) + (rho[2,3] * ((Y2FH
> +
> (delta[2,2] * tau[2,2] * W2FH)) / population[3,1]))))) - ((mew +
> gamma[3,2])
> * W1MM);
> 		dW2MH.dt = (VMH * (partner[3,1] * beta[2,1] * ((rho[3,1] * ((Y2FL +
> (delta[2,2] * tau[2,2] * W2FL)) / population[1,1])) + (rho[2,2] * ((Y1FM
> +
> (delta[2,2] * tau[2,2] * W2FM)) / population[2,1])) + (rho[3,3] * ((Y2FH
> +
> (delta[2,2] * tau[2,2] * W2MH)) / population[3,1]))))) - ((mew +
> gamma[3,2])
> * W1MH);
> 		return(list(c(dXFL.dt, dXFM.dt, dXFH.dt, dXML.dt, dXMM.dt, dXMH.dt,
> dY1FL.dt, dY1FM.dt, dY1FH.dt, dY1ML.dt, dY1MM.dt, dY1MH.dt, dY2FL.dt,
> dY2FM.dt, dY2FH.dt, dY2ML.dt, dY2MM.dt, dY2MH.dt, dZFL.dt, dZFM.dt,
> dZFH.dt,
> dZML.dt, dZMM.dt, dZMH.dt, dVFL.dt, dVFM.dt, dVFH.dt, dVML.dt, dVMM.dt,
> dVMH.dt, dW1FL.dt, dW1FM.dt, dW1FH.dt, dW1ML.dt, dW1MM.dt, dW1MH.dt,
> dW2FL.dt, dW2FM.dt, dW2FH.dt, dW2ML.dt, dW2MM.dt, dW2MH.dt)));
> 	})
> }
> 
> #giving the parameters
> 
> mew = 1/15 																													#proportion of individuals
> entering
> or exiting the sexually active group at a time
> total = 60020																												#total population of sexually
> active
> phi = 0.9 																													#Proportion of individuals who are
> successfully vaccinated
> sigma = 1/10																												#loss of vaccination status
> gamma = matrix(data=c(0.66, 0.66, 0.66, 0.66, 0.66, 0.66, 0.66, 0.66),
> ncol=2, nrow=4)																#Duration of infectiousness
> omega = matrix(data=c(0.82, 0.15, 0.03))
> #proportion
> of those in each sexual activity level
> population = matrix(data=c(omega[1,1]*0.5*total, omega[2,1]*0.5*total,
> omega[3,1]*0.5*total, omega[1,1]*0.5*total, omega[2,1]*0.5*total,
> omega[3,1]*0.5*total), ncol=2, nrow=3)	#population in each sexual
> activity
> level and gender group
> partner = matrix(data=c(1.4, 3, 9))																									#average
> number
> of partners per year by risk group
> beta = matrix(data=c(0.8, 0.7))																									#Chance of
> infection
> from partner given 1 sexual encounter
> rho = matrix(data=c(1, 1, 1, 1, 1, 1, 1, 1, 1), ncol=3, nrow=3)
> #preference for selecting sexual partner by activity group
> tau = matrix(data=c(1, 1, 1, 1), ncol=2, nrow=2)
> delta = matrix(data=c(1, 1, 1, 1), ncol=2, nrow=2)
> #reduction of infection from a vaccinated infected individual
> 
> XFL0 = 1850						#Initial number of females in the low risk group
> unvaccinated susceptibles
> XFM0 = 340						#Initial number of susceptible unvaccinated females
> medium
> risk
> XFH0 = 70						#Initial number of susceptible unvaccinated females high
> risk
> XML0 = 1850						#Initial number of susceptible unvaccinated males low
> risk
> XMM0 = 340						#Initial number of susceptible unvaccinated males medium
> risk
> XMH0 = 70						#Initial number of susceptible unvaccinated males high
> risk
> Y1FL0 = 590						#Initial number of infected unvaccinated females low
> risk
> Y1FM0 = 100						#Initial number of infected unvaccinated females medium
> risk
> Y1FH0 = 20						#Initial number of infected unvaccinated females high
> risk
> Y1ML0 = 590						#Initial number of infected unvaccinated males low risk
> Y1MM0 = 100						#Initial number of infected unvaccinated males medium
> risk
> Y1MH0 = 20						#Initial number of infected unvaccinated males high risk
> Y2FL0 = 590						#Initial number of infected unvaccinated females low
> risk
> Y2FM0 = 100						#Initial number of infected unvaccinated females medium
> risk
> Y2FH0 = 20						#Initial number of infected unvaccinated females high
> risk
> Y2ML0 = 590						#Initial number of infected unvaccinated males low risk
> Y2MM0 = 100						#Initial number of infected unvaccinated males medium
> risk
> Y2MH0 = 20						#Initial number of infected unvaccinated males high risk
> ZFL0 = 0						#Initial number of immune females low risk
> ZFM0 = 0	 					#Initial number of immune females medium risk
> ZFH0 = 0 						#Initial number of immune females high risk
> ZML0 = 0 						#Initial number of immune males low risk
> ZMM0 = 0 						#Initial number of immune males medium risk
> ZMH0 = 0 						#Initial number of immune males high risk
> VFL0 = 21070					#Initial number of susceptible vaccinated females low
> risk
> VFM0 = 3850						#Initial number of susceptible vaccinated females medium
> risk
> VFH0 = 770						#Initial number of susceptible vaccinated females high
> risk
> VML0 = 21070					#Initial number of susceptible vaccinated males low risk
> VMM0 = 3850						#Initial number of susceptible vaccinated males medium
> risk
> VMH0 = 770						#Initial number of susceptible vaccinated males high risk
> W1FL0 = 1110					#Initial number of infected vaccinated females low risk
> W1FM0 = 200						#Initial number of infected vaccinated females medium
> risk
> W1FH0 = 40						#Initial number of infected vaccinated females high risk
> W1ML0 = 1110					#Initial number of infected vaccinated males low risk
> W1MM0 = 200						#Initial number of infected vaccinated males medium risk
> W1MH0 = 40						#Initial number of infected vaccinated males high risk
> W2FL0 = 1110					#Initial number of infected vaccinated females low risk
> W2FM0 = 200						#Initial number of infected vaccinated females medium
> risk
> W2FH0 = 40						#Initial number of infected vaccinated females high risk
> W2ML0 = 1110					#Initial number of infected vaccinated males low risk
> W2MM0 = 200						#Initial number of infected vaccinated males medium risk
> W2MH0 = 40						#Initial number of infected vaccinated males high risk
> 
> p = list(mew=mew, total=total, phi=phi, sigma=sigma, gamma=gamma,
> omega=omega, population=population, partner=partner, beta=beta, rho=rho,
> tau=tau, delta=delta)
> 
> y0 = c(XFL0, XFM0, XFH0, XML0, XMM0, XMH0, Y1FL0, Y1FM0, Y1FH0, Y1ML0,
> Y1MM0, Y1MH0, Y2FL0, Y2FM0, Y2FH0, Y2ML0, Y2MM0, Y2MH0, ZFL0, ZFM0, ZFH0,
> ZML0, ZMM0, ZMH0, VFL0, VFM0, VFH0, VML0, VMM0, VMH0, W1FL0, W1FM0,
> W1FH0,
> W1ML0, W1MM0, W1MH0, W2FL0, W2FM0, W2FH0, W2ML0, W2MM0, W2MH0)
> 
> #Running the ode integrator
> 
> steps= 10;
> t = seq(from=0, to=100, by=1);
> 
> Is anyone able to help?
> 
> 
> 
> 
> 
> --
> View this message in context:
> http://r.789695.n4.nabble.com/DeSolver-giving-NA-as-output-but-running-fully-tp4706497.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

____________________________________________________________
FREE 3D EARTH SCREENSAVER - Watch the Earth right on your desktop!



More information about the R-help mailing list