topic How does JMP compute AICc? in Discussions

How does JMP compute AICc?

nikles — Fri, 09 Jun 2023 22:11:17 GMT

Hi. I'd like to know how JMP computes AICc in the Distribution platform, or in the Fit Censored function? Research on wiki suggests it should be this format:

AICc = 2k - 2*Ln(maxlikelihood) + (2k^2 + 2k)/(n - k - 1)

where n = # datapoint, k = # fit parameters.

However, in this topic another 2 different formulas are provided.

AICc2 = n*Ln(sse/n) + 2*k + (2*k^2 + 2*k)/(n - k - 1);

where sse is the sum of the squares of there residuals between the data and the fitted distribution

And the other formula in the same topic:

AICc3 = n*Ln(sse/n) + (n*(n+k))/(n-k-2);

I've found none of these formulas to agree with each other. To check their consistency, I implemented this script:

Names Default to Here(1);

dt = New Table("Test", 
	Add Rows(100),
	New Column("Random Normal Data", Numeric,
		Formula(Random Normal(0, 1))
	)
);
dt << Run Formulas;
Wait(0.001);
:Random Normal Data << Suppress Eval;		//So data doesn't change between computations.

//Compute AICc from Distribution platform
dist = dt << Distribution(Column("Random Normal Data"));
dist << Fit Distribution(Normal(Density Curve(1)));				//AICc = +285

//Independently compute AICc using method described here: https://en.wikipedia.org/wiki/Akaike_information_criterion
x_vec = :Random Normal Data << Get Values;
mu = Mean(x_vec);
sig = Stddev(x_vec);
maxlikelihood = Normal Distribution(mu, mu, sig);		//e.g. 0.5
k=2;
n=100;
aicc1 = 2*k - 2*Ln(maxlikelihood) + (2*k^2 + 2*k)/(n - k - 1);	//AICc = +5.51
Show(aicc1);

//Independently compute AICc using method described in JMP Forum: https://community.jmp.com/t5/Discussions/how-do-you-get-an-r2-and-AICc-value-for-a-custom-nonlinear/m-p/187375/highlight/true#M40622
ord = Rank(x_vec);
xsort = x_vec[ord];
p_vec = (1::n)`/n;
sse = SSQ(p_vec - Normal Distribution(xsort, mu, sig));
aicc2 = n*Ln(sse/n) + 2*k + (2*k^2 + 2*k)/(n - k - 1);			//AICc = -763
Show(aicc2);

//Independently compute AICc using Other method described in JMP Forum: https://community.jmp.com/t5/Discussions/how-do-you-get-an-r2-and-AICc-value-for-a-custom-nonlinear/m-p/9715
aicc3 = n*Ln(sse/n) + (n*(n+k))/(n-k-2);						//AICc = -661
Show(aicc3);

My results:

Distribution Platform AICc = +285

Formula from wiki AICc = +5.51

Fomula1 in post AICc = -763

Formula2 in post AICc = -661

Since the initial data is generated from a random function, your results may vary somewhat. Nonetheless, it should be apparent that none of these formulas agrees with any other. Most importantly, they do not agree with the value given in the distribution platform.

So, does anyone know how the Distribution platform computes AICc? I need this because I've developed some of my own distributions and I'd like to compare these to those existing in JMP.

FTR, using JMP Pro 15.2.1 on a Mac with Catalina OS.

Thanks.

Re: How does JMP compute AICc?

txnelson — Mon, 19 Apr 2021 03:55:52 GMT

Here is the direct link to the Help on AIC which I got to by selecting the ? tool and then clicking on the AIC table

The AIC Help

Re: How does JMP compute AICc?

Mark_Bailey — Tue, 20 Apr 2021 17:33:17 GMT

The calculation of likelihood and SSE are incorrect, so the second, third, and fourth AICc are also incorrect. The third and fourth formulas use different approximations of the penalty, so will not agree.

Please see these current formulas, which agree with the article in Wikipedia.

Names Default To Here( 1 );

n = 100;

Random Reset( 27513 );
y = J( n, 1, Random Normal() );

dt = New Table( "Test AICc Formula", New Column( "Random Normal Data", Numeric, Continuous, Values( y ) ) );

// Compute AICc with Distribution platform
dist = dt << Distribution( Y( :Random Normal Data ) );
dist << Fit Distribution( Normal( Density Curve( 1 ) ) );
distr = dist << Report;
aiccFromDistribution = distr["Parameter Estimates"][Number Col Box( 4 )] << Get( 2 );
Show( aiccFromDistribution );
dist << Close Window;

// Compute AICc with Fit Least Squares platform
fit = dt << Fit Model(
	Y( :Random Normal Data ),
	Effects,
	Personality( "Standard Least Squares" ),
	Emphasis( "Minimal Report" ),
	Run(
		:Random Normal Data << {Summary of Fit( 1 ), Analysis of Variance( 1 ),
		Parameter Estimates( 1 ), Lack of Fit( 0 ), Scaled Estimates( 0 ),
		Plot Actual by Predicted( 0 ), Plot Regression( 0 ),
		Plot Residual by Predicted( 0 ), Plot Studentized Residuals( 0 ),
		Plot Effect Leverage( 0 ), Plot Residual by Normal Quantiles( 0 ),
		Box Cox Y Transformation( 0 ), AICc( 1 )}
	)
);
fitr = fit << Report;
aiccFromOLS = fitr["Summary of Fit"][NumberColBox(2)] << Get( 1 );
sse = fitr["Analysis of Variance"][NumberColBox(2)] << Get( 2 );
Show( aiccFromOLS );
fit << Close Window;

// Compute AICc using method described in JMP Help (https://www.jmp.com/support/help/en/16.0/#page/jmp/likelihood-aicc-and-bic.shtml#ww293087)
loglikelihood = Summation( i = 1, n, Log( Normal Density( y[i] ) ) );
k = 2;
aiccFromHelp = -2*loglikelihood + 2*k + (2*k*(k+1)/(n-(k+1)));
Show( aiccFromHelp );

//Independently compute AICc using method described in JMP Forum: https://community.jmp.com/t5/Discussions/how-do-you-get-an-r2-and-AICc-value-for-a-custom-nonlinear/m-p/187375/highlight/true#M40622
sse = Sum( (y - Mean( y )) ^ 2 );
aiccFromForum = n*Log( sse/n ) + 2*k + (2*k*(k+1)/(n-(k+1))) + n*Log( 2*Pi() ) + n;
Show( aiccFromForum );