Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Help for AIC and REML score #114

Open
iandanilevicz opened this issue Apr 15, 2024 · 1 comment
Open

Help for AIC and REML score #114

iandanilevicz opened this issue Apr 15, 2024 · 1 comment

Comments

@iandanilevicz
Copy link

iandanilevicz commented Apr 15, 2024

Hi all,

I have a problem with the AIC when I make a pffr model. I fitted a model like this:

pffr_adj_fit1 = pffr(MIMS ~ age + gender + BMI, data = nhanes_ave[1:500,])

and I found the AIC value like this:

pffr_adj_fit1$aic

However, how is this value calculated? Because I found an astronomical value that doesn't coincide with the formula presented by Krivobokova (2007) or Crainiceanu (2024).

n=500

n*log(pffr_adj_fit1$deviance) + 2*sum(pffr_adj_fit1$edf) # Krivobokova equation 4

log(pffr_adj_fit1$deviance) + 2*sum(pffr_adj_fit1$edf)/n # Crainiceanu page 42

I have a second issue regarding how the REML score is calculated. I can find it, but it doesn't match any proxy that I have.

pffr_adj_fit1$gcv.ubre

References:
Krivobokova, T., & Kauermann, G. (2007). “A Note on Penalized Spline Smoothing with Correlated Errors.” Journal of the American Statistical Association, 102(480), 1328–1337, doi: 10.1198/016214507000000978

Crainiceanu, C.M. et al. (2024). "Functional Data Analysis with R". CRC Press, doi: 10.1201/9781003278726

@jaimelopez98
Copy link

Hi,
I tried to achieve the values of AIC and GCV with a simple pffr model using the notation introduced in the help (gam {mgcv}) without success. My code is:

n = 2684
mod1 = pffr(y ~x, method = "GCV.Cp"
mod1$aic #143591443
mod1$deviance / n+ 2*mod1$scale.estimated*sum(mod1$edf)/ n- mod1$scale.estimated # 252615346
mod1$gcv.ubre #65786.27
(n * mod1$deviance) /((n-sum(mod1$edf))^2) #268332499

To check my reference, it is available in the R:
?mgcv::gam

PS. Unfortunately, I didn't solve the problem.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants