如何计算条件期望 Weibull 模型？
Posted
技术标签:
【中文标题】如何计算条件期望 Weibull 模型？【英文标题】：How to calculate the conditional expectation Weibull model? 【发布时间】：2021-03-21 02:16:27 【问题描述】：我想计算 Weibull 模型的条件期望。具体来说，我想估计一个客户在其总任期中随机时刻（时间 = t）的剩余任期。
为此，我计算了每个客户（当前活跃或不活跃）的总任期，并根据每个客户的随机时刻计算他/她当时的任期。
下面的示例是我尝试的快照。我使用 2 个变量 STED 和 TemporalTenure 来预测因变量 tenure，其状态为 0 = 活动或 1 = 非活动。我使用survival包获取生存对象（km_surv）。
df = structure(list(ID = c(16008, 21736, 18851, 20387, 30749,
42159), STED = c(2,
5, 1, 3, 2, 2), TemporalTenure = c(84, 98, 255, 392, 108, 278
), tenure = c(152, 166, 273, 460, 160, 289), status = c(0, 0,
1, 0, 1, 1)), row.names = c(NA,
6L), class = "data.frame")
km_surv <- Surv(time = df$tenure, event = df$status)
df <- data.frame(y = km_surv, df[,!(names(df) %in% c("tenure","status", "ID"))])
weibull_fit <- psm(y ~. , dist="weibull", data = df)
quantsurv <- Quantile(weibull_fit, df)
lp <- predict(weibull_fit, df, type="lp")
print(quantsurv(0.5, lp))
这些估计的输出太高了。我认为这是由包含 TemporalTenure 引起的，但我不知道 psm 包是如何计算这个的，以及是否有其他包可以估计客户 i 的剩余任期时间t。
如何获得以客户已经处于活动状态的时间为条件的预测任期（随机时间：TemporalTenure），其中依赖的任期可以是仍然活跃的客户或非活跃的客户？
编辑
为了澄清，每当我添加时间条件变量，例如：TemporalTenure、收到的付款数量和直到时间 t 的投诉数量时，预测的生命周期在许多情况下都会爆炸。因此，我怀疑psm 不是正确的方法。类似的问题被问到here，但由于同样的原因，给出的解决方案不起作用。
低于已经导致问题的稍大的数据集。
df = structure(list(ID= c(16008, 21736, 18851, 20387, 30749,
42159, 34108, 47511, 47917, 61116, 66600, 131380, 112668, 90799,
113615, 147562, 166247, 191603, 169698, 1020841, 1004077, 1026953,
1125673, 1129788, 22457, 1147883, 1163870, 1220268, 2004623,
1233924, 2009026, 2026688, 2031284, 2042982, 2046137, 2043214,
2033631, 2034252, 2068467, 2070284, 2070697, 2084859, 2090567,
2087133, 2087685, 2095100, 2095720, 2100482, 2105150, 2109353,
28852, 29040, 29592, 29191, 31172, 2126369, 2114207, 2111947,
2102678, 237687, 1093221, 2111607, 2031732, 2105275, 2020226,
1146777, 1028487, 1030165, 1098033, 1142093, 1186763, 2005605,
2007182, 2021092, 2027676, 2027525, 2070471, 2070621, 2072706,
2081862, 2085084, 2085353, 2094429, 2096216, 2109774, 2114526,
2115510, 2117329, 2122045, 2119764, 2122522, 2123080, 2128547,
2130005, 30025, 24166, 61529, 94568, 70809, 159214), STED = c(2,
5, 1, 3, 2, 2, 3, 1, 2, 2, 2, 2, 2, 1, 2, 2, 4, 1, 4, 3, 2, 4,
1, 1, 2, 1, 4, 1, 1, 1, 2, 4, 2, 5, 4, 1, 4, 2, 5, 3, 2, 1, 4,
2, 1, 5, 3, 1, 1, 5, 2, 2, 2, 2, 3, 4, 3, 5, 1, 1, 5, 2, 5, 1,
3, 5, 3, 1, 1, 1, 2, 2, 2, 2, 1, 2, 1, 3, 5, 2, 2, 1, 2, 1, 2,
3, 1, 1, 3, 5, 1, 2, 2, 2, 2, 1, 2, 1, 3, 1), TemporalTenure = c(84,
98, 255, 392, 108, 278, 120, 67, 209, 95, 224, 198, 204, 216,
204, 190, 36, 160, 184, 95, 140, 256, 142, 216, 56, 79, 194,
172, 155, 158, 78, 24, 140, 87, 134, 111, 15, 126, 41, 116, 66,
60, 0, 118, 22, 116, 110, 52, 66, 0, 325, 323, 53, 191, 60, 7,
45, 73, 42, 161, 30, 17, 30, 12, 87, 85, 251, 120, 7, 6, 38,
119, 156, 54, 11, 141, 50, 25, 33, 3, 48, 58, 13, 113, 25, 18,
23, 2, 102, 5, 90, 0, 101, 83, 44, 125, 226, 213, 216, 186),
tenure = c(152, 166, 273, 460, 160, 289, 188, 72, 233, 163,
266, 266, 216, 232, 247, 258, 65, 228, 252, 99, 208, 324,
201, 284, 124, 84, 262, 180, 223, 226, 146, 92, 208, 155,
202, 179, 80, 185, 64, 184, 120, 65, 6, 186, 45, 120, 170,
96, 123, 12, 393, 391, 64, 259, 73, 42, 69, 141, 47, 229,
37, 19, 37, 17, 155, 99, 319, 188, 75, 11, 49, 187, 180,
55, 52, 209, 115, 93, 88, 6, 53, 126, 31, 123, 26, 26, 24,
9, 114, 6, 111, 4, 168, 84, 112, 193, 294, 278, 284, 210),
status = c(0, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 1,
0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 1,
0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 1, 1, 0,
1, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 0, 0, 0, 1, 0, 1), TotalValue = c(2579.35, 2472.85,
581.19, 2579.35, 2472.85, 0, 1829.18, 0, 936.79, 2098.2,
850.47, 2579.35, 463.68, 463.68, 2171.31, 3043.03, 561.16,
3043.03, 3043.03, -68.06, 2098.2, 2504.4, 1536.67, 2719.7,
3043.03, 109.91, 2579.35, 265.57, 3560.34, 2266.95, 3123.16,
3544.4, 1379.19, 2288.35, 2472.85, 2560.48, 1414.45, 3741.49,
202.2, 2856.23, 1457.75, 313.68, 191.32, 2266.95, 661.01,
0, 2050.81, 298.76, 1605.44, 373.86, 3043.03, 2579.35, 448.63,
3043.03, 463.68, 977.28, 818.06, 2620.06, 0, 3235.8, 280.99,
0, 0, 194.04, 3212.75, -23.22, 1833.46, 1829.18, 2786.7,
0, 0, 3250.38, 936.79, 0, 1045.21, 3043.03, 1988.36, 2472.85,
1197.94, 0, 313.68, 3212.75, 1419.33, 531.14, 0, 96.28, 0,
142.92, 174.79, 0, 936.79, 156.19, 2472.85, 463.68, 3520.69,
2579.35, 3328.87, 2567.88, 3043.03, 1081.14)), row.names = c(NA,
100L), class = "data.frame")
【问题讨论】：
所以赏金只是个玩笑，你从来没有打算兑现这个提议？对我来说没那么重要，但它可能会影响人们如何评估你的声誉。 对不起，我希望仍然有不同的解决方案，因为您的解决方案不会改变我完整数据集的结果。我不得不同意，我很惊讶你的解决方案有不同的结果，但是当我在完整的数据集上运行它时，我仍然获得了很高的任期。因此，感觉 psm 解决方案不是为有条件的任期而设计的。 疯狂估计的原因是在协变量中包含 ID 列。 （而且您还没有更正上面代码中可能存在的错误。）我还展示了其他一些不恰当的做法，可以遵循更好的做法。如果没有清楚地了解您在“完整数据集”中看到的内容，我看不出进一步的努力会有多成果。 我很抱歉。我现在已经发布了编辑中的前 100 行，其中包含一个额外的变量“TotalValue”。如果我运行您的解决方案，我会在一种情况下得到 4500 的结果。这大约是最大实际任期的 10 倍。当我像这样添加更多变量时，这些任期变得更加与众不同。psm 函数构造一个带有可能很长的尾巴的参数函数。它很可能会预测很长时间。我还认为从名为TemporalTenure 的预测器预测tenure 存在潜在问题。在 status==0 组中，您提供的数据几乎是一条完美的直线。
【参考方案1】：
这就是我所做的：1) 添加了库调用以加载 pkg:rms，删除了将 Surv 对象放在数据框列中的尝试，3) 在公式中构建 Surv 对象，因为 Therneau 期望公式要构建，并从最可能不属于的协变量中删除 ID。
library(survival); library(rms)
#km_surv <- Surv(time = df$tenure, event = df$status)
#df <- data.frame(y = km_surv, df[,!(names(df) %in% c("tenure","status"))])
weibull_fit <- psm(Surv(time = tenure, event = status) ~TemporalTenure +STED , dist="weibull", data = df)
quantsurv <- Quantile(weibull_fit, df)
lp <- predict(weibull_fit, df, type="lp")
结果#
print(quantsurv(0.5, lp))
1 2 3 4 5 6
151.4129 176.0490 268.4644 466.8266 164.8640 301.2630
【讨论】：
以上是关于如何计算条件期望 Weibull 模型？的主要内容，如果未能解决你的问题，请参考以下文章