# NMaximize runs endlessly - help!

I am doing a maximum likelihood estimation of a parameter (called beta). My likelihood function code is below, and below that I put all the other code which is mostly initialization and the estimator function. It takes in a value of beta and spits out the probability that the data (simdata) would be observed given the value of beta that was fed into the function.

When I run: NMaximize[{lfxn[x], 0 < x && x < 3}, {x}] it runs endlessly.

If I run: Plot[Log[lfxn[x]], {x, 0, 3}] it plots a parabolic curve with a very clear maximum.

Why can't I use NMaximize in the way I currently have to find the max?

## CODE

lfxn[x_] := (
beta = x;
likelihood = 1;
tmatrix =
Table[0, {n + 1}, {n + 1}]; (*Initialize transition matrix*)
Do[tmatrix[[i, i]] = (1 - (Gain[i - 1] + Loss[i - 1])), {i, 2, n}];
Do[tmatrix[[i, i - 1]] = (Loss[i - 1]), {i, 2, n}];
Do[tmatrix[[i, i + 1]] = (Gain[i - 1]), {i, 2, n}];
tmatrix[[1, 1]] = 1;
tmatrix[[n + 1, n]] = Loss[n];
tmatrix[[n + 1, n + 1]] = (1 - Loss[n]);

probvector =
Dot[p, MatrixPower[tmatrix, time]]; (*Probability distribution of states @ time t = V0       dot Tmatrix^(t)*)


(For each value of # of infected in the simulation, determine the likelihood and aggregate accross all simulations for a combined likelihood)

For[i = 1, i <= Length[simdata], i++,
likelihood *= probvector[[simdata[[i]]]]];
likelihood
)


Rest of the code:

time = 1000;
\$RecursionLimit = time + 100;
infected = 2;
dtt = .01;
n = 100;
b = 0;
gamma = .25;
beta = .56;
Gain[i_] := (beta*i*(n - i)*dtt)/n; (*Formula from Brauer*)
Loss[i_] := (b + gamma)*i*dtt;(*Formula from Brauer*)
tmatrix = Table[0, {n + 1}, {n + 1}]; (*Initialize transition matrix*)
probmatrix = Table[0, {time + 1}, {n + 1}];
p = Table[0, {n + 1}]; (*Initial state vector*)
p[[infected + 1]] = 1;
probmatrix[[1, 3]] = 1;
Do[tmatrix[[i, i]] = (1 - (Gain[i - 1] + Loss[i - 1])), {i, 2, n}];
Do[tmatrix[[i, i - 1]] = (Loss[i - 1]), {i, 2, n}];
Do[tmatrix[[i, i + 1]] = (Gain[i - 1]), {i, 2, n}];
tmatrix[[1, 1]] = 1;
tmatrix[[n + 1, n]] = Loss[n];
tmatrix[[n + 1, n + 1]] = (1 - Loss[n]);
dtmcsim = DiscreteMarkovProcess[p, tmatrix];
simnum = 24;
simdata = RandomVariate[dtmcsim[time], simnum];
likelihoodarray = {};
estmin = .01;
estmax = 3;
step = .01;


(Create an array of elements of the form {a,b} where a = parameter \ estimate, b = likelihood of observed data given a)

For[k = estmin, k < estmax, k += step,
likelihoodarray = Append[likelihoodarray, {k, lfxn[k]}]
];

-