Re: Precision/FullSimplify
- To: mathgroup at smc.vnet.net
- Subject: [mg56280] Re: Precision/FullSimplify
- From: "Jens-Peer Kuska" <kuska at informatik.uni-leipzig.de>
- Date: Thu, 21 Apr 2005 05:35:59 -0400 (EDT)
- Organization: Uni Leipzig
- References: <d458t7$i0o$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Hi, what many Chop[] do ?? Regards Jens "leigh pascoe" <leigh at cephb.fr> schrieb im Newsbeitrag news:d458t7$i0o$1 at smc.vnet.net... > Gentle beings, > > I hope the method of cutting and pasting from > Mathematica shown below is > acceptable. If not I would appreciate > instructions on how to do it in a > more acceptable fashion. > Consider the likelihood equations for two > parameters x and y, based on > experimental data (a,b,c,d,e,f,g) that are > observed numbers (integers). > The likelihood and derivatives are defined in Ma > by: > > L:=(b+d+f)*Log[x]+(e+g)*Log[y]-(a+b)*Log[1+x]-(f+g)*Log[x+y]-(c+d+e)*Log[1+2*x+y] > D[L,x] > D[L,y] > > Yielding the likelihood equations: > > eq1=(b+d+f)/x-(a+b)/(1+x)-2*(c+d+e)/(1+2*x+y)-(f+g)/(x+y)\[Equal]0 > eq2=(e+g)/y-(c+d+e)/(1+2*x+y)-(f+g)/(x+y)\[Equal]0 > > which I would like to solve for the maximum > likelihood estimates of the > parameters and then evaluate for different data > sets. Firstly define the > information matrix: > > m11:=-D[D[L,x],x] > m12:=-D[D[L,x],y] > m21:=-D[D[L,y],x] > m22:=-D[D[L,y],y] > m:={{m11,m12},{m21,m22}}//MatrixForm > m > > which can be inverted when the equations are > solved to give the > variances of the estimates. To solve the > equations I tried > > Sols=Solve[{eq1,eq2},{x,y}] > > Which gives me the answers I want, however the > expressions are quite > large. By inspection I find that the last > solution is the one that > corresponds to the maximum likelihood. However I > run into problems > evaluating the resulting expression. If I try to > evaluate the last > solution numerically I get a complex number, I > assume due to rounding > errors. e.g. > > N[Sols] > /.{a->50,b->50,c\[Rule]2,d\[Rule]4,e->100,f\[Rule]2,g->100} > > produces > > \!\({{y \[Rule] \(-5.336646740150529`\) + > 1.2546145287986817`*^-14\ > \[ImaginaryI], > x \[Rule] > \(\(3.8558775093812896`\)\(\[InvisibleSpace]\)\) > - > 2.220446049250313`*^-16\ > \[ImaginaryI]}, {y \[Rule] \ > \(-1.067199413695643`\) - > 3.6306168234410265`*^-14\ \[ImaginaryI], > x \[Rule] \(-0.4135698170735974`\) - > 6.661338147750939`*^-16\ > \[ImaginaryI]}, {y \[Rule] \(\(50.`\)\(\ > \[InvisibleSpace]\)\) + > 1.1823875212257917`*^-14\ \[ImaginaryI], > x \[Rule] \(\(1.`\)\(\[InvisibleSpace]\)\) > + > 6.661338147750939`*^-16\ > \[ImaginaryI]}}\) > > The last solution is the correct one but without > the addition of the > small complex part. I then tried > > N[Sols, $MaxPrecision] > /.{a->50,b->50,c\[Rule]2,d\[Rule]4,e->100,f\[Rule]2, > g->100} > > Which gives me presumably exact expressions, but > it is still not clear > that they are real. So I tried > > FullSimplify[ N[Sols, $MaxPrecision] > /.{a->50,b->50,c\[Rule]2,d\[Rule]4,e->100,f\[Rule]2, > g->100}] > N[%] > > {{y\[Rule]-5.33665,x\[Rule]3.85588},{y\[Rule]-1.0672, > x\[Rule]-0.41357},{y\[Rule]50.,x\[Rule]1.}} > > and the last solution appears to be the one I > want, but it takes a long > time. Now let's try this with a different data > set. > > FullSimplify[ > N[Sols, $MaxPrecision] > /.{a\[Rule]42,b\[Rule]90,c\[Rule]5,d\[Rule]29, > e\[Rule]49,f\[Rule]6,g\[Rule]17}] > > This one seems to hang forever. In any case > there must be a better way > of evaluating the solutions for different data > sets. Should I be telling > Mathematica to look for real solutions, using > NSolve, substituting in > eq1 and eq2 before solving?? It is preferable to > generate all the > solutions before choosing the one that gives the > maximum likelihood. Any > suggestions, observations, help will be greatly > appreciated. > > Leigh > (Inexperienced User) >