在2015年的二月21日，我的妻子已经33天没有来月经了，她怀孕了，这真是天大的好消息！通常月经的周期是大约一个月，如果你们夫妇打算怀孕，那么月经没来或许是一个好消息。但是33天，这还无法确定这是一个消失的月经周期，或许只是来晚了，那么它是否真的是一个好消息？

为了能获得结论我建立了一个简单的贝叶斯模型，基于这个模型，可以根据你当前距离上一次经期的天数、你历史经期的起点数据来计算在当前经期周期中你怀孕的可能性。在此篇文章中我将阐述我所使用的数据、先验思想、模型假设以及如何使用重点抽样法获取数据并用R语言运算出结果。在最后，我将解释为什么模型的运算结果最终并不重要。另外，我将附上简便的脚本以供读者自行计算.

非常幸运的是，在2014年的下半年间我的妻子一直在记录她经期起始日期，否则我只能以仅拥有小量数据而告终。总体上我们拥有8个经期的起始日期数据，但是我采用的数据不是日期而是相邻经期起始日间相隔的天数。 已经有33天。

period_onset <- as.Date(c("2014-07-02", "2014-08-02", "2014-08-29", "2014-09-25",
                          "2014-10-24", "2014-11-20", "2014-12-22", "2015-01-19"))
days_between_periods <- as.numeric(diff(period_onset))

所以日期发生得相对规律，以28天为一个周期循环。最后一次月经开始日期是在1月19日，所以在2月21日，距离最后一次经期发生日。

我要建立一个涵盖生理周期的模型，包括受孕期和不受孕期，这显然需要做大量的简化。我做了一些总体假设如下：

接下来是我所做的具体假设：

基本的假设就是这样了。但是为了使其更加实际，需要考虑使用一个似然函数，一个给定了参数和一些数据、计算在给定参数下数据的概率，通常而言是一个与概率成正比例的数值——似然值。因为这个似然值可能极小所以我需要对其取对数，从而避免引起数值问题。当用R语言设计似然函数时，总体上的模式如下：

calc_log_like <- function(days_since_last_period, days_between_periods,
                          mean_period, sd_period, next_period,
                          is_fertile, is_pregnant) {
  n_non_pregnant_periods <- length(days_between_periods)
  log_like <- 0
  if(n_non_pregnant_periods > 0) {
    log_like <- log_like + sum( dnorm(days_between_periods, mean_period, sd_period, log = TRUE) )
  }
  log_like <- log_like + log( (1 - 0.19 * is_fertile)^n_non_pregnant_periods )
  if(!is_pregnant && next_period < days_since_last_period) {
    log_like <- -Inf
  }
  log_like
}

这里数据有标量days_since_last_period以及向量days_between_periods，而其他的参数将会被被估计出来。使用这个函数，我能从任意一个数据+参数的组合中得出对数似然函数值。但是，到这里我只完成了建模的一半工作，我还需要先验信息！

为了完善这个模型，我需要所有参数的先验信息。换言之，我需要明确在获取数据之前这个模型包含了哪些信息。具体上，我需要实验开始前mean_period, sd_period, is_fertile, and is_pregnant的初始值。（虽然next_period也是一个参数，我不需要给出一个它的确切初始值，因为它的分布完全由mean_period 和sd_period确定。另外，我还需要找到在一个周期内能受孕的可能值（上文中我设定为0.19）。这里我使用了模糊、主观的数据吗？不！我到生育文献中去寻找了更加有信息价值的依据！

对于days_between_periods的分布，其参数为mean_period和sd_period。这里我使用了来自文章The normal variabilities of the menstrual cycle Cole et al, 2009 中的估计值，该文测量了184个年龄来自18-36岁的女性的经期规律。相邻经期间天数的总平均值为27.7天。每一个参与实验者的标准差的平均值为2.4。总体样本的间隔天数的标准差为1.6。给定了这些估计值以后我令mean_period服从（27.7,2.4）的正态分布，令sd_period服从均值为1.6，标准差为2.05的半正态分布。如下：

对于参数is_fertile a以及参数is_pregnant我考虑了受精频率作为先验。想要确定可育的夫妻的比例几乎是不可能的事情，因为这里对于不育有各种不同的定义。 Van Geloven et al. (2013)做了一个小范围的文献回顾然后得出结论所有夫妻中有2%至5%的人被认为是不孕的。因为曾看到高达10%的情况，我决定取该范围的上限。设定初始数据100%-5%=95%的夫妻是可孕的。

is_pregnant 是 0 1变量表示这对夫妻在最近的一轮周期中是否将要（或者说已经）受孕。在这里我使用的先验值是在一个周期内成功受孕的概率。当这对夫妇没有生育能力时这个概率值显然为0.0，但是积极地尝试、可育的夫妇在一个周期内成功受孕的比例有多大呢？不幸的是我并没有找到明确说明这一数据的文献，但是我找到了比较接近的参照依据。在Increased Infertility With Age in Men and Women Dunson et al. (2004) 一书的第53页，给出了在12个月中一直尝试受孕但是没有怀上的夫妻的比例，同时该数据也提供了女性不同年龄段的数据。