I have co-opted some R code that estimates parameters for a Gompertz-Makeham mortality model. The data to which the model is fit is interval censored, with counts of individuals in each age-range category. For example:

Code:

`naga <- structure(c(15,20,35,50,20,35,50,Inf,46,310,188,101),.Dim=as.integer(c(4,3)),.Dimnames=list(NULL,c("col1","col2","col3")))`

Code:

```
GM.naga <- function(x,deaths=naga)
{
a2=x[1]
a3=x[2]
b3=x[3]
shift<-15
nrow<-NROW(deaths)
S.t<-function(t)
{
return(exp(-a2*(t-shift)+a3/b3*(1-exp(b3*(t-shift)))))
}
d<-S.t(deaths[1:nrow,1])-S.t(deaths[1:nrow,2])
obs<-deaths[,3]
lnlk<-as.numeric(crossprod(obs,log(d)))
return(lnlk)
}
optim(c(0.001, 0.01, 0.1),GM.naga,control=list(fnscale=-1))
```

Code:

```
surv_GM.naga <- function (t)
{
x=c(5.227131e-09, 2.208408e-02, 4.577932e-02)
a2<-x[1]
a3<-x[2]
b3<-x[3]
shift<-15
S.t<-exp(-a2*(t-shift)+a3/b3*(1-exp(b3*(t-shift))))
return<-S.t
print(S.t)
}
surv_GM.naga(15:70)
```

--Trey