2014-10-24 6 views
6

Ho fatto un semplice esperimento di un GLM in statsmodels, ed è stato perplesso scoprire perché i risultati GLM non contenevano attributi R^2?Perché i modelli di statistiche GLM non hanno R^2 nei risultati?

Mi sembra che qui ci sia qualcosa di molto semplice sul perché GLM non abbia un calcolo R^2, e modi in cui posso calcolarlo da solo.

Grazie!

In [1]: import pandas as np 

In [2]: import pandas as pd 

In [3]: import numpy as np 

In [4]: import statsmodels.api as sm 

In [5]: data = pd.DataFrame({'col1':np.arange(10),'col2':np.arange(
KeyboardInterrupt 

In [5]: x = np.arange(0,10,0.5) 

In [6]: 

In [6]: y = np.zeros(len(x)) 

In [7]: y[0] = 0 

In [8]: for i in range(1,len(x)): 
    ...:   y[i] = 0.5*x[i] + 2.5*y[i-1] + 10*np.random.rand() 
    ...:  

In [9]: print y 
[ 0.00000000e+00 9.35177024e-01 8.18487881e+00 2.95126464e+01 
    8.08584645e+01 2.11423251e+02 5.38685230e+02 1.35653420e+03 
    3.39564225e+03 8.49234338e+03 2.12377817e+04 5.31015961e+04 
    1.32764789e+05 3.31924691e+05 8.29818265e+05 2.07455796e+06 
    5.18640343e+06 1.29660216e+07 3.24150658e+07 8.10376747e+07] 

In [10]: X = pd.DataFrame({'x1':x[1:],'y-Lag1':y[:-1]}) 


In [11]: m1 = sm.GLM(y[1:],X).fit() 

In [12]: m1.summary() 
Out[12]: 
<class 'statsmodels.iolib.summary.Summary'> 
""" 
       Generalized Linear Model Regression Results     
============================================================================== 
Dep. Variable:      y No. Observations:     19 
Model:       GLM Df Residuals:      17 
Model Family:    Gaussian Df Model:       1 
Link Function:    identity Scale:     12.9022715725 
Method:       IRLS Log-Likelihood:    -50.199 
Date:    Thu, 23 Oct 2014 Deviance:      219.34 
Time:      13:44:22 Pearson chi2:      219. 
No. Iterations:      3           
============================================================================== 
       coef std err   t  P>|t|  [95.0% Conf. Int.] 
------------------------------------------------------------------------------ 
x1    1.5746  0.175  8.999  0.000   1.232  1.918 
y-Lag1   2.5000 1.23e-07 2.03e+07  0.000   2.500  2.500 
============================================================================== 
""" 

In [13]: m1. 
m1.aic     m1.llf     m1.remove_data 
m1.bic     m1.load     m1.resid_anscombe 
m1.bse     m1.model     m1.resid_deviance 
m1.conf_int    m1.mu      m1.resid_pearson 
m1.cov_params    m1.nobs     m1.resid_response 
m1.deviance    m1.norm_resid    m1.resid_working 
m1.df_model    m1.normalized_cov_params m1.save 
m1.df_resid    m1.null     m1.scale 
m1.f_test     m1.null_deviance   m1.summary 
m1.family     m1.params     m1.summary2 
m1.fit_history   m1.pearson_chi2   m1.t_test 
m1.fittedvalues   m1.pinv_wexog    m1.tvalues 
m1.initialize    m1.predict     
m1.k_constant    m1.pvalues   
+0

Curioso! prova a stampare il risultato di m1.rsquared()? – user3378649

+0

Controlla anche questa discussione: http://stackoverflow.com/questions/21785049/statsmodels-poisson-glm-different-than-r – user3378649

+0

Non sei sicuro del motivo per cui è stato implementato in questo modo, ma la risposta SO sopra e alcuni wikipedia-ing mi ha permesso di calcolare a mano la R^2: sst_val = sum (map (lambda x: np.power (x, 2), y-np.mean (y))) sse_val = sum (map (lambda x: np. power (x, 2), m1.resid_response)) r2 = 1.0 - sse_val/sst_val –

risposta

1

Non è sicuro circa il motivo per cui la sua implementato in questo modo, ma il SO risposta di cui sopra e un po 'wikipedia-ing mi ha permesso di calcolare la R^2 a mano facilmente:

sst_val = sum(map(lambda x: np.power(x,2),y-np.mean(y))) 
sse_val = sum(map(lambda x: np.power(x,2),m1.resid_response)) 
r2 = 1.0 - sse_val/sst_val 
4

Per GLM con errori gaussiana e il collegamento dell'identità, R^2 ha senso (se il modello ha una costante), ma non ha senso come una generale bontà della misura adatta per GLM. È possibile presentare una richiesta di miglioramento (o, meglio ancora, una richiesta di pull) per includere alcune migliori e più generali informazioni sulle statistiche di adattamento nei risultati GLM.

È possibile leggere ulteriori informazioni al riguardo nel contesto di R here.

Problemi correlati