Go back to Flu Mortality

library(deSolve)

dyn.load(“SIR.so”)

data1 = c(106, 121, 136, 156, 183, 248, 279, 276, 239, 185, 160, 187, 150, 193, 153, 153, 138, 143, 119)

t1 = c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18)

model_SIRn = function(t, theta) {

SIR = function(Time, State, Pars) {
with(as.list(c(State, Pars)), {

	DecreaseS = -1 * beta * S * I / N
	DecreaseI = -1 * alpha * I

	dS = DecreaseS
	dI = -1 * DecreaseS + DecreaseI

	return(list(c(dS, dI)))
	})
}	

with(as.list(theta), {
	
state1 = c(S = N - I1, I = I1)	
pars1 = c(beta = beta1, alpha = alpha1, N = N)
	
out1 = ode(y = state1, times = t, func = SIR, parms = pars1, 	method = 'ode45')
	
	return(-diff(out1[,"S"]))
})

}

model_SIRc = function(t, theta) {

SIR = function(Time, State, Pars) 
{
	dS = double(1)
	dI = double(1)
	
         out = .C("SIR",
                 Time = as.integer(Time),
                 State = as.double(State),
                 Pars = as.double(Pars),
                 dS = as.double(dS),
                 dI = as.double(dI))
                 
         return(list(c(out$dS, out$dI)))
}

with(as.list(theta), {
	
state1 = c(S = N - I1, I = I1)	
pars1 = c(beta = beta1, alpha = alpha1, N = N)
	
out1 = ode(y = state1, times = t, func = SIR, parms = pars1, 	method = 'ode45')

	return(-diff(out1[,"S"]))
})

}

theta = c(beta1 = 0.3, alpha1 = 0.2, N = (sum(data1)*1.166), I1 = 1)

norm = model_SIRn(c(0:19), theta)

print(“normal one”) print(norm)

c = model_SIRc(c(0:19), theta)

print(“c one”) print©