如何在scipy.stats.gamma.fit中获取拟合参数的误差估计？

1 回答

• 7

  edit 警告：以下说明了问题中示例后面 GenericLikelihoodModel 的用法 . 然而，在伽马分布的情况下，位置参数改变了分布的支持，该支持被最大似然估计的一般假设排除 . 更标准的用法将修复支持，floc = 0，因此它只是一个双参数分布 . 在这种情况下，标准MLE理论适用 .

  Statsmodels有一个用于最大似然估计的通用类， GenericLikelihoodModel . 它不是直接针对这种情况设计的，但可以在一些帮助下使用（定义属性并提供start_params） .

  import numpy as np
  from statsmodels.base.model import GenericLikelihoodModel
  
  from scipy.stats import gamma
  shape = 12; loc = 0.71; scale = 0.0166
  data = gamma.rvs(shape, loc=loc, scale=scale, size=1000)
  params = gamma.fit(data) # params close to but not the same as (shape, loc, scale) 
  # HOW TO ESTIMATE/GET ERRORS FOR EACH PARAM?
  
  print(params)
  print('\n')
  
  
  class Gamma(GenericLikelihoodModel):
  
      nparams = 3
  
      def loglike(self, params):
          return gamma.logpdf(self.endog, *params).sum()
  
  
  res = Gamma(data).fit(start_params=params)
  res.df_model = len(params)
  res.df_resid = len(data) - len(params)
  print(res.summary())
  

  这将打印以下内容

  (10.31888758604304, 0.71645502437403186, 0.018447479022445423)
  
  
  Optimization terminated successfully.
           Current function value: -1.439996
           Iterations: 69
           Function evaluations: 119
                                  Gamma Results                                 
  ==============================================================================
  Dep. Variable:                      y   Log-Likelihood:                 1440.0
  Model:                          Gamma   AIC:                            -2872.
  Method:            Maximum Likelihood   BIC:                            -2852.
  Date:                Sun, 12 Jul 2015                                         
  Time:                        04:00:05                                         
  No. Observations:                1000                                         
  Df Residuals:                     997                                         
  Df Model:                           3                                         
  ==============================================================================
                   coef    std err          z      P>|z|      [95.0% Conf. Int.]
  ------------------------------------------------------------------------------
  par0          10.3187      2.242      4.603      0.000         5.925    14.712
  par1           0.7165      0.019     37.957      0.000         0.679     0.753
  par2           0.0184      0.002      8.183      0.000         0.014     0.023
  ==============================================================================
  

  然后，基于最大似然估计的其他结果也是可用的，例如，可以通过指定限制矩阵或通过使用具有相等性的字符串表达式来执行第一参数为10的z测试：

  >>> res.t_test(([1, 0, 0], [10]))
  <class 'statsmodels.stats.contrast.ContrastResults'>
                               Test for Constraints                             
  ==============================================================================
                   coef    std err          z      P>|z|      [95.0% Conf. Int.]
  ------------------------------------------------------------------------------
  c0            10.3187      2.242      0.142      0.887         5.925    14.712
  ==============================================================================
  
  >>> res.t_test('par0=10')
  <class 'statsmodels.stats.contrast.ContrastResults'>
                               Test for Constraints                             
  ==============================================================================
                   coef    std err          z      P>|z|      [95.0% Conf. Int.]
  ------------------------------------------------------------------------------
  c0            10.3187      2.242      0.142      0.887         5.925    14.712
  ==============================================================================
  

  回复于 2024-04-28T11:38:44+08:00