?p=18782

本文我们讨论了期望寿命的计算。人口统计模型的起点是死亡率表。但是，这种假设有偏差，因为它假设生活条件不会得到改善。为了正确处理问题，我们使用了更完整的数据，其中死亡人数根据x岁而定，还包括日期t。


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      DE=read.table(
      "DE.txt",skip = 
      3,header=
      TRUE)
     
    
   
    
     
    
    
     
      EXPS=read.table(
      "EXPS.txt",skip = 
      3,header=
      TRUE)

我们用 Dx，t表示死亡人数，Ex，t表示暴露人数。因此，对于在日期t上x岁的某人，在该年死亡的概率为 qx，t = Dx，t / Ex，t。这些数据存储在矩阵中进行可视化，存储在数据库中进行回归。

必须进行一些修改以避免出现零值的问题，因为（i）我们求出比率（ii）然后我们对数化）。我们可以可视化为x和t的函数。

persp(log(QF))

或


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      persp3d(
      ages,annees,log(
      QH),col=
      "light blue")

为了模拟qx，t的演化，我们可以从Lee＆Carter（1992）的模型中获得启发，该模型假设log （qx，t）= Ax + Bx⋅Kt。A =（A0，A1，⋯，A110）在某种程度上是log（qx，t）。K =（K1816，K1817，⋯，K2015）使我们了解生活条件的改善，一年内死亡的可能性降低。这些改善不是均匀的，因此我们使用B =（B0，B1，⋯，B110）来使改善取决于l '年龄。

为了估计参数A，B和K，我们尝试使用二项式模型。B（Ex，t，qx，t），这是人寿保险的基本模型。这里Dx，t〜B（Ex，t，exp [ Ax + Bx⋅Kt]）。

另一个线索是使用小数定律，即如果概率低（一年中的死亡概率就是这种情况），则二项式定律可以近似由泊松分布。我们在这里用到了Poisson回归，其解释变量为年龄x，年t和暴露量为偏移变量。唯一的问题是它不是线性回归。我们这里有非线性模型，因为E [ Dx，t] =（exp[log（Ex，t）+ Ax + Bx⋅Kt]）。


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      gnm( 
      DH ~ offset(
      log(
      EH)  + as.factor(
      age) +
     
    
   
    
     
    
    
     
      Multas.factor(
      age,as.factor(
      annee),
     
    
   
    
     
    
    
     
      family = poisson(
      link=
      "log")

我们有估计系数A ^，B ^和K ^。


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      Ax=reg
      $coefficients[
      2:
      111]
     
    
   
    
     
    
    
     
      Bx=reg
      $coefficients[
      112:
      222]
     
    
   
    
     
    
    
     
      Kt=reg
      $coefficients[
      223:length(reg
      $coefficients)]

我们可以表示三组系数。首先 A ^表示平均变化，

plot(ages[-1],Ax)

我们还可以用 K ^来绘制时间。

同样，该模型不可被识别。简而言之，改善没有任何意义。我们可以表示-K ^，它的优点是描述了生活条件的改善。最后，让我们作图-B ^

困难在于，为了预测期望寿命，我们需要针对t的大值（尚未观察到）计算qt，x。例如，某人可能想知道q50,2020（对于1970年出生的人）。我们要使用q50,2020 = exp（A ^ 50 + B ^ 50 K ^ 2020）。问题是K ^ 2020不属于估计数量K ^。

这个想法是Lee＆Carter（1992）的初衷，我们可以尝试指数模型或线性模型（在1950年以后的原始K ^序列上）


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      lm(log(Kt[idx])~ann[idx])
     
    
   
    
     
    
    
     
      futur=
      2016:
      2125
     
    
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      lm(Kt[idx]~ann[idx])
     
    
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      points(futur,pr,col=
      "blue")

然后，我们可以根据过去的数据建立一系列预测，q ^ x，t = exp [A ^ x + B ^ x K ^ t]，以及未来数据q〜x，t = exp [A ^ x + B ^ x K〜t]。

我们保留过去的数据，这里是1880年死亡的概率


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      plot(BASE
      $x[BASE
      $t==1880],BASE
      $pred[BASE
      $t==1880],
     
    
   
    
     
    
    
     
      log=
      "y")

同样，我们在未来（此处为2050年）使用这两种模型


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      BASE2$Qpred
      1=exp(cste+BASE
      2$Ax+BASE
      2$Bx*BASE
      2$Kt
      1)
     
    
   
    
     
    
    
      
     
    
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      plot(BASE
      2$x[BASE
      2$t==
      2050],BASE
      2$Qpred
      1[BASE
      2$t==
     
    
   
    
     
    
    
     
      2050],log=
      "y")

用于指数预测

对于线性预测，对1968年出生的人，我们有第二年死亡的概率


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      if(sbase$t[i]<= 
      2015)
     
    
   
    
     
    
    
     
      {v
      q[i]=BASE[ BASE$x==sbase$x[i]) &  BASE$t==sbase$t[i]),
      "Qpred"] 
     
    
   
    
     
    
    
     
      if(sbase$t[i] <
      2015) 
     
    
   
    
     
    
    
     
      {v
      q[i]=BASE2[(BASE2$x==sbase$x[i]) & (BASE2$t==sbase$t[i]),
      "Qpred2"]

左边是我们模型估算值，右边是预测值。

要计算出生时的期望寿命，我们使用以下代码


  
   
    
     
    
    
     
      sum(cumprod(
      exp(-v
      q[1:110])))
     
    
   
    
     
    
    
     
      [
      1] 
      77.62047

然后，我们可以做函数可视化这种期望寿命的演变


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      vP = cumprod(exp(-(sbase$vq[
      1:
      110])))
     
    
   
    
     
    
    
     
      sum(vP)}


  
   
    
     
    
    
      
     
    
   
    
     
    
    
     
      ANN =
      1930:
      2010
     
    
   
    
     
    
    
     
      plot(ANN ,E
      2)