How to fit a poisson glm in R when the parameters are not separated I am trying to fit a poisson GLM to count data.  
Where the count $N_r$ for time $t_r$ is assumed to follow $\mathrm{ Poisson} (\lambda (t_r - t_0)^{-\alpha} e^{-\beta t_r})$.  
Here $r$ is discrete and we observe $N_r$ and $t_r$ for $r$ from 1 to 30. The parameters that I want to fit are $\lambda, \alpha, \beta, t_0 $. Intuitively I am trying to see in what form the counts $N_r$ decrease as time $t_r$ increases.  
In R, I can fit the model as follows in the absence of $t_0$ (this is not what I wanted to do, but what i have been able to do) as follows.
$$ model1 = glm(N_r \sim t_r + \log(t_r), family = poisson()) $$ 
Since,  $\log(E(N_r))$ is a linear combination $\log(\lambda)$, $t_r$, and $\log(t_r)$.  
The results of model1 is reasonable.  
However, I am having difficulty when trying to include the $t_0$ term in my model equation. I would like to fit a model that would estimate the parameters $\lambda, \alpha, \beta, t_0 $. How can i do this in R? It is confusing for me because it appears that,
$$ \log(E(N_r)) = \log(\lambda) - \beta t_r - \alpha \log(t_r - t_0)$$ 
So $\alpha$ and $t_0$ are not separated. (Perhaps I have to use a different function?) It would be great to get some insight on how to fit this model.
 A: Unless you have just the right combination of $\alpha$, $\beta$, and values of $x_i$, it may be impossible to pin down the value of $t_0$.  The problem is that these curves can be  difficult to distinguish.
To obtain a fit, consider searching over a reasonable set of values of $t_0$.  For each value the problem is a standard Poisson Generalized Linear Model, and so can be fit routinely.  Extract a measure of the fit, such as deviance (or equivalent measures AIC and BIC).  A "plausible" value of $t_0$ is any one whose fit isn't too much greater than the value with smallest deviance.  A good cutoff threshold to use for the difference is twice an upper percentile of the chi-squared distribution with 1 degree of freedom.  The range of resulting values of $t_0$ is a confidence interval for $t_0$.
By choosing suitable units of measurement for $t_r$ (which I will simply call $t$), we may suppose the sample values of $t$ range between $0$ to $1$.  This makes it easier to compare models and gives $t_0$ an inherent meaning.  For instance, $t_0\approx -1$ is fairly large because it's comparable to the full range of the $t$. We will limit our search to negative values of $t_0$, because the model makes no sense when $t_0$ equals or exceeds the smallest $t$ in the dataset.
Here is an example generated from the model with $\log(\lambda)=0$, $\alpha=2$, $\beta=1/3$, and $t.0 = -1/10$, with $30$ $t$ values equally spaced throughout its range, as described in the question.

I searched a range of $t.0$ values geometrically spaced from $-0.01$ to $-3$. Here is the plot of how the deviance varies:

Points within half the 95th percentile of $\chi^2(1)$ (equal to $1.92$) are denoted "significant."  The corresponding values of $t_0$ range from $-0.022$ to $-0.37$.  Such a narrow range of estimates is rare: in many cases the deviance will scarcely change, no matter what value is given to $t_0$.
The next figure helps show why estimates of $t_0$ can range so widely.  For each one of the "significant" values of $t_0$ in the search, I plot the curve corresponding to the GLM fit of $\alpha$, $\beta$, and $\log(\lambda)$.  You can see the curves hardly vary:

They vary the most at low values of $t$.  This is valuable information.  It indicates that if you had a choice concerning which values of $t$ to measure, you should concentrate many of them at small values. As an example, I created a dataset in which $t$ was sampled at $6$ geometrically spaced values from $0.007$ to $1$.  Each one of these values was measured $5$ times independently, producing another dataset of $30$ values. Here are the plausible fits along with the data:

There is almost no visible difference among them.  The range of plausible values of $t_0$ is noticeably narrower, demonstrating the improvement in estimation that is afforded by a better choice of values of $t$ to sample.
To obtain confidence intervals for $\alpha$, $\beta$, and $\log(\lambda)$, extract the confidence intervals from the GLM model fit with $t_0$ set to the one with smallest deviance.
The fits manage to remain the same by rebalancing the effects of the other parameters as $t_0$ is varied.  Using this first dataset again, here is how the estimates vary with $t_0$ (as $t_0$ ranges through its $95\%$ confidence interval--the "significant" values shown in the second figure).

(The true value of $t_0$ is marked with a vertical line.)
Even if you restrict the models to non-negative values of $\alpha$ and $\beta$, there's still quite a range of $(\log(\lambda),\alpha,\beta,t_0)$ parameter combinations that produce almost the same fit.

For those interested in pursuing this further, here is the R code that produced the figures.
library(ggplot2)
#
# Specify the data-generation model.
#
lambda.log <- 0
t.0 <- -0.1
alpha <- 2
beta <- 1/3
size <- 0.95 # Percentile of chi-square(1) for deviance cutoff
# n <- 6
# x <- rep(exp(seq(-5, 0, length.out=n)), 5)
x <- seq(0, 1, length.out=30)
f <- function(x, theta) {
  exp(theta["lambda.log"] - log(x - theta["t.0"]) * theta["alpha"] - theta["beta"] * x)
}
theta <- c(lambda.log=lambda.log, alpha=alpha, beta=beta, t.0=t.0)
#
# Generate data.
#
set.seed(17)
X <- data.frame(x=x, y.hat=f(x, theta))
X$y <- with(X, rpois(length(x), y.hat))
#
# Generate a closely-spaced set of regressors for the purpose of plotting
# fits, etc.
3
x.0 <- seq(min(x), max(x), length.out=201)
X.0 <- data.frame(x=x.0, y.hat=f(x.0, theta))
#
# FIGURE 1: Data with the underlying model.
#
ggplot(X, aes(x,y)) + geom_path(aes(x,y.hat), X.0, size=1) + 
  geom_point(alpha=1/2, shape=21, fill="Black") + #geom_smooth() + 
  xlab("t") + 
  ggtitle("Data with True Curve")
#
# Search a range of values of t.0.
#
u <- -exp(seq(-5, 1, length.out=31))
fits <- lapply(u, function(u) glm(y ~ I(log(x-u)) + x, X, family="poisson"))
#
# Identify which values of t.0 are closest to achieving the best deviance.
#
deviance <- sapply(fits, function(x) x$deviance)
i <- which(deviance <= min(deviance) + qchisq(size, 1)/2)
#
# FIGURE 2: Deviance vs. t.0
#
D <- data.frame(t.0=u, deviance=deviance, 
                Significant=deviance <=  min(deviance) + qchisq(size, 1)/2)
ggplot(D, aes(x=t.0, y=deviance)) + 
  geom_line(size=1) + 
  geom_point(aes(color=Significant, shape=Significant), size=2) + 
  ggtitle("Deviance vs. t.0")
#
# FIGURE 3: Display plausible fits.
#
m <- length(i)
X.1 <- data.frame(t.0=rep(u[i], each=length(x.0)),
                  x = rep(x.0, m), 
                  y.hat=c(sapply(i, function(j) {
                    b <- c(lambda.log=1, alpha=-1, beta=-1) * fits[[j]]$coef
                    f(x.0, c(b, t.0=u[j]))})))
X.1$t.0 <- ordered(signif(X.1$t.0, 1))
ggplot(X.1, aes(x,y.hat)) + 
  geom_line(aes(color=t.0, group=t.0), size=1, alpha=1/2) + 
  geom_point(aes(x,y), X, alpha=1/2, shape=21, fill="Black") +
  xlab("t") + ylab("y") + 
  ggtitle("Plausible Fits For Varying t.0")
#
# FIGURE 4: Plot coefficients vs. t.0.
#
library(data.table)
a <- sapply(fits[i], coef)
A <- as.data.table(t(a))
names(A) <- c("lambda.log", "alpha", "beta")
A$t.0 <- u[i]
A <- rbind(A[, .(t.0, Estimate=lambda.log, Coefficient="log(lambda)")],
           A[, .(t.0, Estimate=-alpha, Coefficient="alpha")],
           A[, .(t.0, Estimate=-beta, Coefficient="beta")])
ggplot(A, aes(t.0, Estimate, group=Coefficient)) + 
  geom_vline(xintercept=t.0) +
  geom_line(aes(color=Coefficient), size=1) + 
  ggtitle("Coefficient Estimates")

A: Even the model code you gave does not do what you seem to say it does. It assumes
$$\log \text{E} N_r = \beta_0 + \beta_1 t_r + \beta_2 \log t_r$$
and estimates $\beta_0$, $\beta_1$ and $\beta_2$, which I believe is not what you had in mind.
You may need to use some software that lets you specify non-linear transformations of parameters or lets you specify a completely user defined likelihood. For the latter there are a lot of options: e.g. PROC NLMIXED in SAS, the rstan R package with it's maximum likelihood capabilities (or perhaps if your parameters are not that well identified you could go Bayesian) and probably a lot of other packages.
