r - 为双变量生存审查数据编写语法以拟合 R 中的 copula 模型

Question

library(Sunclarco)
library(MASS)
library(survival)
library(SPREDA)
library(SurvCorr)
library(doBy)

#Dataset
diabetes=data("diabetes")

data1=subset(diabetes,select=c("LASER","TRT_EYE","AGE_DX","ADULT","TIME1","STATUS1"))
data2=subset(diabetes,select=c("LASER","TRT_EYE","AGE_DX","ADULT","TIME2","STATUS2"))


#Adding variable which identify cluster
data1$CLUSTER<- rep(1,197)
data2$CLUSTER<- rep(2,197)

#Renaming the variable so that that we hve uniformity in the common items in the data
names(data1)[5] <- "TIME" 
names(data1)[6] <- "STATUS"
names(data2)[5] <- "TIME"
names(data2)[6] <- "STATUS"

#merge the files
Total_data=rbind(data1,data2)

# Re arranging the database 
diabete_full=orderBy(~LASER+TRT_EYE+AGE_DX,data=Total_data)
diabete_full

#using Sunclarco package for Clayton a nd Gumbel

Clayton_1step <- SunclarcoModel(data=diabete_full,time="TIME",status="STATUS",
                          clusters="CLUSTER",covariates=c("LASER","TRT_EYE","ADULT"),
                          stage=1,copula="Clayton",marginal="Weibull")
summary(Clayton_1step)

#              Estimates    StandardErrors

#lambda       0.01072631    0.005818201
#rho          0.79887565    0.058942208
#theta        0.10224445    0.090585891
#beta_LASER   0.16780224    0.157652947
#beta_TRT_EYE 0.24580489    0.162333369
#beta_ADULT   0.09324001    0.158931463

#              Estimate      StandardError
#Kendall's Tau 0.04863585    0.04099436



Clayton_2step <- SunclarcoModel(data=diabete_full,time="TIME",status="STATUS",
                                clusters="CLUSTER",covariates=c("LASER","TRT_EYE","ADULT"),
                                stage=2,copula="Clayton",marginal="Weibull")
summary(Clayton_1step)

#             Estimates    StandardErrors
#lambda       0.01131751    0.003140733
#rho          0.79947406    0.012428824
#beta_LASER   0.14244235    0.041845100
#beta_TRT_EYE 0.27246433    0.298184235
#beta_ADULT   0.06151645    0.253617142
#theta        0.18393973    0.151048024

#                Estimate     StandardError
#Kendall's Tau   0.08422381   0.06333791



Gumbel_1step <- SunclarcoModel(data=diabete_full,time="TIME",status="STATUS",
                                clusters="CLUSTER",covariates=c("LASER","TRT_EYE","ADULT"),
                                stage=1,copula="GH",marginal="Weibull")

#             Estimates     StandardErrors
#lambda       0.01794495     0.01594843
#rho          0.70636113     0.10313853
#theta        0.87030690     0.11085344
#beta_LASER   0.15191936     0.14187943
#beta_TRT_EYE 0.21469814     0.14736381
#beta_ADULT   0.08284557     0.14214373

#                Estimate      StandardError
#Kendall's Tau   0.1296931     0.1108534



Gumbel_2step <- SunclarcoModel(data=diabete_full,time="TIME",status="STATUS",
                               clusters="CLUSTER",covariates=c("LASER","TRT_EYE","ADULT"),
                               stage=2,copula="GH",marginal="Weibull")

需要在 R 中为不同的 copula 类拟合 copula 模型，特别是 Gaussian、FGM、Pluckett 和可能的 Frank（如果我还有时间）。正在使用的数据是 R 中通过 Survival 和 Survcorr 包提供的糖尿病数据。

这是我的论文正在研究，它是一项探索性研究，以了解 copula 类如何服务于不同的目的，因为它们导致在相同的结果上产生不同的结果。我在 Rstudio 中找到了一个包 Sunclarco，我能够适应 Clayton 和 Gumbel copula 类，但它还不适用于其他类。

面临的挑战是，由于我已经审查了必须纳入可能性估计的数据，因此我编写语法变得更加困难，因为我没有强大的编程背景。此外，我必须合并编程中存在的协变量，看看它们是否存在对关联的影响。然而，他向我的发起人提出了关于如何处理这个难题的语法写作的见解，如下所示

• ******首先，忘记似然函数。我们只使用对数似然函数。这样，您不需要对每个观测值取贡献的乘积，但可以对不同观测值取对数贡献的总和。

• 接下来，由于我们有一个平衡的设计，我们可以使用常规的数据帧结构，其中每个集群在数据帧中只有一行。生命周期、指标和所有协变量等不同变量是该数据框中的列。

• 由于双变量设置，只有 4 种可能的方式可以对对数似然函数做出贡献：均未删失、均删失、第一个未删失和第二个删失，或第一个删失和第二个未删失。好吧，要创建对数似然函数，您需要在数据框中创建一个新变量，在其中根据夫妻中的哪个人被审查，对数似然做出正确的贡献。当你取这个变量的总和时，你就有了对数似然函数的值。

• 由于此函数取决于参数，因此您可以使用任何优化器，例如 optim 或 nlm 来获得最佳值。这里要小心， optim 和 nlm 寻找函数的最小值，而不是最大值。这很容易解决，因为函数 -f 的最小值与函数 f 的最大值相同。

• 由于每个copula 函数都有不同的导数表达式，现在应该可以得到似然函数。******

我仍在努力寻找一种方法，因为每个 copula 类的每个可能性都会发生变化，因为生成器函数对于相应的 copula 也是唯一的，因为它需要在估计期间进行调整。最后，我应该对 copula 估计的一个和两个步骤进行分析，因为我将使用它来比较结果。

如果有人可以帮助我解决这个问题，那么我将永远感激不尽。即使只有一个 copula 类，例如 Gaussian，那么我会根据请求帮助的那个来计算其余部分，因为我尝试了所有方法，但我仍然没有任何东西可以出现，现在我觉得时间已经不多了我自己。

标签： rstatistics