tweedie 的有限混合

Posted 2023-03-12

【中文标题】tweedie 的有限混合

【英文标题】：Finite mixture of tweedie

【发布时间】：2015-06-05 09:00:34

【问题描述】：

我正在尝试估计 tweedie（或复合 Poisson-gamma）分布的有限混合。我搜索了我能想到的所有资源，但没有找到任何关于如何做到这一点的资源。

我目前正在尝试使用 R 中的 flexmix 包编写不同的 M-step 驱动程序，如第 12-14 页的 flexmix 小插图中所述。这是我的代码，它依赖于 cplm 包：

tweedieClust <- function(formula = .~.,offset = NULL)
require(tweedie)
require(cplm)
require(plyr)
require(dplyr)
retval <- new("FLXMC", weighted = TRUE, formula = formula, dist = "tweedie",
              name = "Compound Poisson Clustering")

retval@defineComponent <- expression (
    predict <- function(x, ...) 
        pr <- mu
    
    logLik <- function(x, y, ...)
        dtweedie(y, xi = p, mu = mu, phi = phi) %>%
             log
    
    new("FLXcomponent",
        parameters=list(coef=coef),
        logLik=logLik, predict=predict,
        df=df)
)
retval@fit <- function (x, y, w, component) 
    fit <- cpglm(formula = y ~ x, link = "log", weights=w, offset=offset)
    with(list(coef = coef(fit), df = ncol(x),mu = fit$fitted.values,
              p = fit$p, phi = fit$phi),
         eval(retval@defineComponent))

retval

但是，这会导致以下错误：

dtweedie(y, xi = p, mu = mu, phi = phi) 中的错误： 不相容数组的二元运算

有没有人做过或见过 tweedie 分布的有限混合？您能否指出我使用 flexmix 或其他方式完成此任务的正确方向？

【问题讨论】：

您找到解决方案了吗？

【参考方案1】：

问题出在权重部分，如果你删除它，它会起作用：

tweedieClust <- function(formula = .~.,offset = NULL)
  require(tweedie)
  require(statmod)
  require(cplm)
  require(plyr)
  require(dplyr)
  retval <- new("FLXMC", weighted = F, formula = formula, dist = "tweedie",
            name = "Compound Poisson Clustering")

  retval@defineComponent <- expression (
    predict <- function(x, ...) 
      pr <- mu
    
    logLik <- function(x, y, ...)
      dtweedie(y, xi = p, mu = mu, phi = phi) %>%
        log
    
    new("FLXcomponent",
        parameters=list(mu=mu,xi=p,phi=phi),
        logLik=logLik, predict=predict,
        df=df)
  )
  retval@fit <- function (x, y, w, component) 
    fit <- cpglm(formula = End~.,data=dmft, link = "log")
    with(list(df = ncol(x), mu = fit$fitted.values,
              p = fit$p, phi = fit$phi),
         eval(retval@defineComponent))
  
  retval

示例：

library(flexmix)
data("dmft", package = "flexmix")
m1 <- flexmix(End ~ .,data=dmft, k = 4, model = tweedieClust())