How to calculate AIC(Akaike's information criterion) for restricted cubic spline

[Nobuaki Michihata]

I would like to calculate the AICs (Akaike's information criterion) to find out which one is more appropriate for different restricted cubic spline with different numbers of knots.
Does anyone know how to calculate the AICs?

Thank you for your cooperation.

The followings are sample commands.

capture ssc install postrcspline

sysuse uslifeexp, clear

mkspline2 ys = year, cubic nknots(3)
des ys*
reg le ys*

adjustrcspline, ///
addplot(scatter le year, msymbol(Oh)) ///
ytitle("life expectancy (knots3)")
graph export rcsknots3.tif, as(tif) replace

capture drop ys*
mkspline2 ys = year, cubic nknots(5)
des ys*
reg le ys*

adjustrcspline, ///
addplot(scatter le year, msymbol(Oh)) ///
ytitle("life expectancy (knots5)")
graph export rcsknots5.tif, as(tif) replace

[Paul Dickman]

Here's one approach. I'm sure there are more efficient ways but hopefully this will give you a start. It also shows how you can add the AIC to your graphs.

The approach is to loop over the desired values of the number of knots. We store the model estimates (estimate store) so we can produce a table of AIC values when we exit the loop. While inside the loop, we write the AIC to a local macro variable so we can, for example, include it with the graph.

sysuse uslifeexp, clear

forvalues i = 3/7 {
capture drop ys*
mkspline2 ys = year, cubic nknots(`i')
quietly regress le ys*
estimates store knots`i'
quietly estat ic
local mAIC: display %4.1f -2*e(ll) + 2*e(rank)    
adjustrcspline, ///
addplot(scatter le year, msymbol(Oh)) title("`i' knots; AIC: `mAIC'") ///
ytitle("life expectancy") name(knots`i',replace)
}

// table of model statistics
estimates table knots*, stats(N r2_a aic bic)

You could also save results to a postfile (help postfile) to give even more flexibility.

There's almost certainly a better place to get the stored AIC without having to use estat ic. I couldn't find it in e() after running -estimates store- and a brief look at the manual didn't help. According to the manual, -estat ic- stores the AIC in the matrix r(S) so you could get it directly from there but you should be able to access it without running -estat ic-. I'll leave it to you to read the manual in detail or maybe someone else will illustrate how to do it.

I rarely use linear regression (I work with categorical and survival outcomes) so don't know all the relevant commands, but I hope this helps.

Thanks for following the forum guidelines and showing a worked example.

[Nobuaki Michihata]

Dear Paul Dickman,

Thank you for your very excellent response.
I'm learning a lot.
I would like to be able to contribute in this forum.
Thank you so much. It was very helpful.