Hola. Me gustaría saber cómo JMP calcula AICc en la plataforma de distribución o en la función Fit Censored. Investigación en wiki sugiere que debería ser este formato:
AICc = 2k - 2 * Ln (probabilidad máxima) + (2k ^ 2 + 2k) / (n - k - 1)
donde n = # punto de datos, k = # parámetros de ajuste.
Sin embargo, en este tema se proporcionan otras 2 fórmulas diferentes.
AICc2 = n * Ln (sse / n) + 2 * k + ( 2 * k ^ 2 + 2 * k) / (n - k - 1 );
donde sse es la suma de los cuadrados de los residuos entre los datos y la distribución ajustada
Y la otra fórmula en el mismo tema:
AICc3 = n * Ln (sse / n) + (n * (n + k)) / (nk- 2 );
No he encontrado ninguna de estas fórmulas que coincida entre sí. Para verificar su consistencia, implementé este 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);Mis resultados:
Plataforma de distribución AICc = +285
Fórmula de wiki AICc = +5.51
Fomula1 en la publicación AICc = -763
Fórmula2 en la publicación AICc = -661
Dado que los datos iniciales se generan a partir de una función aleatoria, sus resultados pueden variar un poco. No obstante, debería ser evidente que ninguna de estas fórmulas concuerda con ninguna otra. Lo más importante es que no están de acuerdo con el valor otorgado en la plataforma de distribución.
Entonces, ¿alguien sabe cómo la plataforma de distribución calcula AICc? Necesito esto porque he desarrollado algunas de mis propias distribuciones y me gustaría compararlas con las que existen en JMP.
FTR, usando JMP Pro 15.2.1 en una Mac con Catalina OS.
Gracias.
Aquí está el enlace directo a la Ayuda sobre AIC al que llegué seleccionando? herramienta y luego haciendo clic en la tabla AIC
La ayuda de AIC
El cálculo de probabilidad y SSE son incorrectos, por lo que el segundo, tercero y cuarto AICc también son incorrectos. Las fórmulas tercera y cuarta usan aproximaciones diferentes de la penalización, por lo que no estarán de acuerdo.
Consulte estas fórmulas actuales, que concuerdan con el artículo de 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 );