Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Trying to solve a MLE with nonlinear constraints that need to be satisfied. However, the given solution do not match the constraints. Am I doing it correctly?

(* Construct Probability Model *)
\[ScriptCapitalD]1 =  NegativeBinomialDistribution[r, \[Alpha]/(\[Alpha] + 1)];
\[ScriptCapitalD]2 =   ParameterMixtureDistribution[BetaBinomialDistribution[a, b, n], n \[Distributed] \[ScriptCapitalD]1];
(* Test for one set of data *)
T = 598; R = 333; A = 78; U = 893; (* X=106; assumed unknown *)
data = {3, 1, 9, 1, 4, 2, 12, 2, 3, 1, 2, 1, 3, 2, 2, 2, 39, 3, 1, 7, 2, 1, 5, 1, 2, 1, 1, 64, 1, 2, 1, 18, 2, 1, 2, 11, 4, 3, 5, 2, 1, 4, 8, 2, 10, 16, 2, 2, 1, 5, 1, 1, 1, 1, 1, 2, 2, 4, 1, 4, 1, 1, 1, 1, 2, 1, 2, 4, 2, 4, 2, 1, 8, 1, 1, 1, 1, 1};
CC = {FullSimplify[Mean[\[ScriptCapitalD]1]] == T/U, Mean[\[ScriptCapitalD]2] == R/U, Probability[z == 0, z \[Distributed] \[ScriptCapitalD]2] == 1 - (A/U)}; (* Constraints *)
CP = {a > 0 && b > 0 && r > 0 && \[Alpha] > 0}; (* Parameter constraints *)
MLE = NMaximize[{LogLikelihood[\[ScriptCapitalD]2, data], Flatten[{CC, CP}]}, {a, b, r, \[Alpha]}]
(* Do they match the constraints? Expressed in terms of known values *)
U*{CC[[1, 1]] , CC[[2, 1]], 1 - CC[[3, 1]]} /. MLE[[2]]
(* X_estimate *)
U*(1 - Probability[y == 0, y \[Distributed] \[ScriptCapitalD]1]) /. MLE[[2]]

When I ran the code received a lot of warnings, including "Simplify::fas: Warning: one or more assumptions evaluated to False." with final output being:

{-338.828, {a -> 4.46969*10^-8, b -> 3.49161*10^-8,  r -> 0.112262, \[Alpha] -> 0.167517}}
{598.446, 335.984, 98.1873}
174.889

However, with trial and error found that the constraints can be satisfied with

{a -> 0.79170676874115343262052276450750130908150706082622,  b -> 0.63003691806728426319650273891989981624150723094542,  r -> 0.045708496520516912401968683076575539420378914371680, \[Alpha] \-> 0.068257002329133115008291066471192704521037320465721}

being one answer. Whether those are the maximum likelihood estimates I am not sure.

share|improve this question

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.