Forum: help

RE: Estimation Parameter Lomax distribution [ Reply ]
By: Arne Henningsen on 2014-05-30 04:48

The warning message "NaNs produced" means that a calculation produced NAs. In your example taking the log or square root of a negative number.

The message "names for target but not for current" means that the object that you specify as argument "target" of function all.equal( ) has names, while the object that you specify as argument "current" does not have names.

The command "all.equal(list(mleHess), coef(mleHess),vcov(mleHess))" produces so many (warning) messages, because you try to compare more than two objects and the compared objects are very different from each other.

RE: Estimation Parameter Lomax distribution [ Reply ]
By: Rian Hidayat on 2014-05-23 00:15

[forum:41073]

Thanks for your reply.Arne.

I drived the gradient and hessians, and :

> data <- c(37506920,16021662,40836717,96892397,15674554,77282747,2635150,389826409,21710782,12808038,13403370,176898918,28375371,27644834,60031467,269584588,35259848,66065117,2232743,1744820,
+ 85942868,6597570,14297220,3826130,18754414,3714655,153941444,22775334,154037993,138826899,157971031,125847655,873180,1763300,3270505,39907212,82304026,34540520,8000000,70658389,5761140,183568978,
+ 87150648,3711400,1088500,15775480,15664513,41260775,99492628,33684222,2864400,20810793,157852579,3926551,207513707,14168105,19449025,31036530,1997187,92499979,37385277,167303755,98294267,30766575,
+ 1251985,5782980,47953884,2751700,95387615,14501130,11508175,39280451,484973165,8511393,91657510,13006543,1575035,45607920,29946778,82065716,22943384,48860523,30592670,14714350,37112564,
+ 44319238,9122050,11055258,65693845,17661581,10424348,12179839,37933080,197711333,10856502,54765307,23483110,97333812,24936761,62507046,52223776,43908853,53847308,2559040,11698815,
+ 29410707,3284558,80479248,5960395,161226418,4115865,28265220,1325205,14397242,40280496,35806655,95270489,3942750,5224756,19008640,31635800,891169,170249292,4457718,188056938,25548340,14153160,223624,
+ 40212318,48851956,105883461,1761331,15093855,13047762,16845752,59383415,166157999,107533281,5498360,123424827,28098035,10306345,1981560,52593805,31860080,53441511,43007545,26910137,39108375,
+ 4815930,161052537,22604038,37078790,55171253,220990,33782809,7491461,10544015,8505884,11205759,24143830,3833060,174697425,19554825,157736242,26552750,14259344,17982270,102757338,12053303,18948790,
+ 6534675,10952615,3844750,108866995,8250603,38274688,35475903,24781428,6005376,8359924,883575,7150048,550900,35875893,61111841,92690529,134803823,7456224,7816480,5419266,14060043,152355304,2669975,33284700,7345310)
> library(maxLik)
> n<-length(data)
> logLikFun=function(pars,x){
+ beta=pars[1]
+ delta=pars[2]
+ n*log(delta)+n*delta*log(beta)-(delta+1)*sum(log(x+beta))
+ }
> mle<-maxLik(logLikFun,start=c(80000000,2), method="NR", x=data)
Warning messages:
1: In log(delta) : NaNs produced
2: In log(delta) : NaNs produced
3: In log(delta) : NaNs produced
4: In log(delta) : NaNs produced
> summary(mle)
--------------------------------------------
Maximum Likelihood estimation
Newton-Raphson maximisation, 7 iterations
Return code 2: successive function values within tolerance limit
Log-Likelihood: -3655.564
2 free parameters
Estimates:
Estimate Std. error t value Pr(> t)
[1,] 8.0000e+07 NA NA NA
[2,] 2.5551e+00 0.0000e+00 Inf < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
--------------------------------------------
Warning messages:
1: In sqrt(diag(vc)) : NaNs produced
2: In sqrt(diag(vc)) : NaNs produced

=======================================
what the meaning "NaNs produced"? and how to solve?
========================================

> library(maxLik)
> logLikGrad <- function( pars,x ) {
+ beta=pars[1]
+ delta=pars[2]
+ N <- length( x )
+ logLikGradValues <- numeric( 2 )
+ logLikGradValues[ 1 ] <- (n*delta/beta)-(delta+1)*sum(1/(x+beta))
+ logLikGradValues[ 2 ] <- (n/delta)+n*log(beta)-sum(log(x+beta))
+ return( logLikGradValues )
+ }
> mleGrad <- maxLik( logLik = logLikFun, grad = logLikGrad,
+ start = c( beta = 80000000, delta = 2 ),x=data )
> all.equal( logLik( mleGrad ), logLik( mle ) )
[1] TRUE
> all.equal( coef( mleGrad ), coef( mle ) )
[1] "names for target but not for current"

=======================================
what the meaning "names for target but not for current"? and how to solve?
========================================

> logLikHess=function(pars,x){
+ beta=pars[1]
+ delta=pars[2]
+ N <- length( x )
+ logLikHessValues=matrix(0,nrow=2, ncol=2)
+ logLikHessValues[1,1]=((-n)*delta/beta^2)+(delta+1)*sum(1/((x+beta)^2))
+ logLikHessValues[1,2]=(n/beta)-sum(1/(x+beta))
+ logLikHessValues[2,1]=(n/beta)-sum(1/(x+beta))
+ logLikHessValues[2,2]=(-n)/delta^2;
+ return(logLikHessValues)
+ }
> mleHess=maxLik(logLikFun, grad=logLikGrad,
+ hess=logLikHess,start = c( beta = 80000000, delta = 2 ),x=data)
> all.equal(list(mleHess), coef(mleHess),vcov(mleHess))
[1] "Modes: list, numeric"
[2] "names for current but not for target"
[3] "Length mismatch: comparison on first 1 components"
[4] "Component 1: Modes: list, numeric"
[5] "Component 1: names for target but not for current"
[6] "Component 1: Attributes: < Modes: list, NULL >"
[7] "Component 1: Attributes: < names for target but not for current >"
[8] "Component 1: Attributes: < Length mismatch: comparison on first 0 components >"
[9] "Component 1: Length mismatch: comparison on first 1 components"
[10] "Component 1: Component 1: Mean relative difference: 21885.45"
Warning messages:
1: In if (check.attributes) attr.all.equal(target, current, ...) :
the condition has length > 1 and only the first element will be used
2: In if (check.attributes) attr.all.equal(target, current, ...) :
the condition has length > 1 and only the first element will be used
3: In if (check.attributes) attr.all.equal(target, current, tolerance = tolerance, :
the condition has length > 1 and only the first element will be used

=======================================
why only the firs element will be used?
========================================

> list(logLik(mleGrad), coef(mleGrad), vcov(mleGrad))
[[1]]
[1] -3655.564

[[2]]
beta delta
8.000000e+07 2.555066e+00

[[3]]
beta delta
beta Inf Inf
delta Inf Inf

Sorry, I have many questions..
and thanks for your Reply.

RE: Estimation Parameter Lomax distribution [ Reply ] By: Arne Henningsen on 2014-05-22 20:54	[forum:41072]
Replace: logLikGrad <- function( pars ) { by: logLikGrad <- function( pars, x ) { and replace: mleGrad <- maxLik( logLik = logLikFun, grad = logLikGrad, start = c( beta = 83551500, delta = 2 ) ) by: mleGrad <- maxLik( logLik = logLikFun, grad = logLikGrad, start = c( beta = 83551500, delta = 2 ), x = data ) ATTENTION: untested code. Best regards, Arne

RE: Estimation Parameter Lomax distribution [ Reply ]
By: Rian Hidayat on 2014-05-22 13:47

[forum:41069]

Thanks for reply

I drived the gradients, and :

> library(maxLik)
> n<-length(data)
> logLikFun=function(pars,x){
+ beta=pars[1]
+ delta=pars[2]
+ n*log(delta)+n*delta*log(beta)-(delta+1)*sum(log(x+beta))
+ }
> res<-maxLik(logLikFun,start=c(83551500,2), method="NR", x=data)
Warning messages:
1: In log(delta) : NaNs produced
2: In log(delta) : NaNs produced
3: In log(delta) : NaNs produced
4: In log(delta) : NaNs produced
5: In log(delta) : NaNs produced
6: In log(delta) : NaNs produced
> summary(res)
--------------------------------------------
Maximum Likelihood estimation
Newton-Raphson maximisation, 7 iterations
Return code 2: successive function values within tolerance limit
Log-Likelihood: -3655.558
2 free parameters
Estimates:
Estimate Std. error t value Pr(> t)
[1,] 8.3552e+07 NA NA NA
[2,] 2.6365e+00 1.7010e-01 15.5 < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
--------------------------------------------
Warning messages:
1: In sqrt(diag(vc)) : NaNs produced
2: In sqrt(diag(vc)) : NaNs produced
>
> logLikGrad <- function( pars ) {
+ beta=pars[1]
+ delta=pars[2]
+ N <- length( x )
+ logLikGradValues <- numeric( 2 )
+ logLikGradValues[ 1 ] <- (n*delta/beta)-(delta+1)*sum(1/(x+beta))
+ logLikGradValues[ 2 ] <- (n/delta)+n*log(beta)-sum(log(x+beta))
+ return( logLikGradValues )
+ }
> mleGrad <- maxLik( logLik = logLikFun, grad = logLikGrad,
+ start = c( beta = 83551500, delta = 2 ) )
Error in x + beta : 'x' is missing

Why x or data is missing?

Thanks

RE: Estimation Parameter Lomax distribution [ Reply ] By: Arne Henningsen on 2014-05-22 09:42	[forum:41068]
You can derive the gradients and let maxLik() use them (see section 3.1 of the paper that I referenced in my previous post) and/or use a different optimisation method.

RE: Estimation Parameter Lomax distribution [ Reply ] By: Rian Hidayat on 2014-05-22 08:26	[forum:41067]
Thank's Arene, For your Replay... I am really sorry, i am not a good programer. so, my log-likehood function cannot obtain to estimation parameter with package maxLik? thanks

RE: Estimation Parameter Lomax distribution [ Reply ]
By: Arne Henningsen on 2014-05-22 05:18

[forum:41062]

Dear Rian

Your log-likelihood function lnl() only returns the likelihood value but not its gradients or its Hessian. Therefore, the gradients and the Hessian are calculated by finite-difference approximations, which can result in small errors. In your case, the upper left element of the Hessian matrix (see "res$hessian") is positive due to approximation errors so that the upper left element of the variance covariance matrix of the estimated coefficients is negative (see "vcov(res)"). As the square root of a negative value is not defined, you cannot obtain the standard error of the first parameter. Please take a look at sections 3.1 and 3.2 of [1].

[1] Henningsen, Arne and Ott Toomet (2011): maxLik: A package for maximum likelihood estimation in R. Computational Statistics 26(3), p. 443-458. http://files.itslearning.com/data/ku/103018/publications/maxlik.pdf

Best regards,
Arne

Estimation Parameter Lomax distribution [ Reply ]
By: Rian Hidayat on 2014-05-22 01:59

[forum:41061]

estimation parameter.docx (54) downloads

Dear Oot & Anne
Thank you so much for making the nice package "maxLik".
I'd liek to find the makimum likelihood estimates of two parameters (in a log-likelihood function) for lomax distribution, but i have a problem when estimation a parameter not exact.
My codes :

> data <- c(37506920,16021662,40836717,96892397,15674554,77282747,2635150,389826409,21710782,12808038,13403370,176898918,28375371,27644834,60031467,269584588,35259848,66065117,2232743,1744820,85942868,6597570,14297220,3826130,18754414,3714655,153941444,22775334,154037993,138826899,157971031,125847655,873180,1763300,3270505,39907212,82304026,34540520,8000000,70658389,5761140,183568978,87150648,3711400,1088500,15775480,15664513,41260775,99492628,33684222,2864400,20810793,157852579,3926551,207513707,14168105,19449025,31036530,1997187,92499979,37385277,167303755,98294267,30766575,1251985,5782980,47953884,2751700,95387615,14501130,11508175,39280451,484973165,8511393,91657510,13006543,1575035,45607920,29946778,82065716,22943384,48860523,30592670,14714350,37112564,44319238,9122050,11055258,65693845,17661581,10424348,12179839,37933080,197711333,10856502,54765307,23483110,97333812,24936761,62507046,52223776,43908853,53847308,2559040,11698815,29410707,3284558,80479248,5960395,161226418,4115865,28265220,1325205,14397242,40280496,35806655,95270489,3942750,5224756,19008640,31635800,891169,170249292,4457718,188056938,25548340,14153160,223624,40212318,48851956,105883461,1761331,15093855,13047762,16845752,59383415,166157999,107533281,5498360,123424827,28098035,10306345,1981560,52593805,31860080,53441511,43007545,26910137,39108375,4815930,161052537,22604038,37078790,55171253,220990,33782809,7491461,10544015,8505884,11205759,24143830,3833060,174697425,19554825,157736242,26552750,14259344,17982270,102757338,12053303,18948790,6534675,10952615,3844750,108866995,8250603,38274688,35475903,24781428,6005376,8359924,883575,7150048,550900,35875893,61111841,92690529,134803823,7456224,7816480,5419266,14060043,152355304,2669975,33284700,7345310)
>
> n<-length(data)
> lnl=function(pars,x){
+ beta=pars[1]
+ delta=pars[2]
+ n*log(delta)+n*delta*log(beta)-(delta+1)*sum(log(x+beta))
+ }
> res<-maxLik(lnl,start=c(83489000,2), method="NR", x=data)
> summary(res)
--------------------------------------------
Maximum Likelihood estimation
Newton-Raphson maximisation, 6 iterations
Return code 2: successive function values within tolerance limit
Log-Likelihood: -3655.558
2 free parameters
Estimates:
Estimate Std. error t value Pr(> t)
[1,] 8.3489e+07 NA NA NA
[2,] 2.6351e+00 0.0000e+00 Inf < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
--------------------------------------------
Warning messages:
1: In sqrt(diag(vc)) : NaNs produced
2: In sqrt(diag(vc)) : NaNs produced

why p-value of parameter beta be NA?
Estimate Std. error t value Pr(> t)
[1,] 8.3489e+07 NA NA NA

thanks