Bayesian Analysis of Single-Molecule Experimental Data

Kou, S. C.; Sunney Xie, X.; Liu, Jun S.

doi:10.1111/j.1467-9876.2005.00509.x

总结

实验技术的最新进展使科学家能够在单分子基础上跟踪生物化学过程，这比传统的集合平均实验提供了更丰富的化学动力学信息，但也带来了许多新的统计挑战。本文首次对用于探测单个脱氧核糖核酸（DNA）发夹分子构象动力学的单分子荧光寿命实验进行了似然统计分析。构象变化最初被视为连续时间的两态马尔可夫链，这是不可观察的，必须从光子发射的变化中推断出来。由于未观测到的分子布朗扩散，该模型更加复杂。除了简单的两态模型外，文献中还提出了一种竞争模型，该模型将DNA发夹的两种状态之间的能量屏障建模为Ornstein–Uhlenbeck过程。我们首先推导了简单两态模型的似然函数，然后推广了处理复杂情况的方法，如未观察到的分子扩散和波动的能量垒。发展了数据增强技术和马尔可夫链蒙特卡罗方法，以从期望的后验分布中进行采样。贝叶斯因子计算和相关参数的后验估计表明，波动障碍模型比简单的两状态模型更符合数据。

贝叶斯因子,布朗扩散,连续时间马尔可夫链,斯过程,能量屏障,可能性,Ornstein–Uhlenbeck工艺,缩放变换更新

1.简介

最近的技术进步使科学家能够对单分子动力学进行观测，这在几十年前是不可想象的（聂和扎尔，1997; 谢和特劳特曼，1998; 韦斯，2000; 塔马拉特等。,2000; 莫纳，2002)-著名物理学家理查德·费曼（Richard Feynman）曾描述说，看到单个原子的图像是一种“宗教体验”。作为对大分子集合进行的传统实验的补充，单分子实验为新的科学发现提供了巨大的潜力和许多优势。首先，可以直接测量分子性质的分布，而不仅仅是系综平均值。其次，单分子实验可以实时跟踪生物化学过程的瞬态中间产物，而这在以前只能通过同步大分子集合的动作来完成。第三，单分子轨道提供了详细的动力学信息，这是传统的系综实验所无法提供的。详细的动态信息对于具有复杂内部结构的复杂生物分子尤其重要（Xie和Lu，1999; 杨等。,2003)因为了解单个分子的动力学对于揭示其生化特性和功能至关重要。

新的进展不仅提供了机遇，也带来了从实验、理论到统计的重大挑战。统计挑战来自单分子轨道的随机性。虽然实验和理论问题吸引了许多研究人员，但单分子研究的统计方面却很少受到关注。我们在这里考虑了一种特殊类型的单分子实验，用于研究脱氧核糖核酸（DNA）发夹动力学，以举例说明几个统计问题，如单分子动力学和实验并发症的建模、感兴趣参数的有效估计、模型判别和统计计算。为了成功应对实验复杂性，我们使用了贝叶斯数据增强方法，该方法提供了比该领域广泛使用的矩型方法更精确的估计。我们在这里采用贝叶斯方法，不仅因为它在处理复杂模型和干扰参数时具有一致性和最佳性，而且因为我们对涉及的各种参数有非常详细的知识（先验知识）。这里开发的一般策略也可以应用于其他单分子实验。

1.1. 脱氧核糖核酸发夹

DNA发夹是一种两端碱基互补的单链核酸结构，因此可以形成分子内配对。由于配对的形成和断裂，DNA发夹有开放和闭合两种状态。在闭合状态下，两端配对在一起，整个结构像发夹，而在开放状态下，分子内配对被打破（阀盖等。,1998; 克里切夫斯基和博内，2002).图1说明了DNA发夹5的闭合和开放状态^′–CCCAA–（T）₂₁–TTGGG–3个^′。

图1

DNA发夹的两种状态：为了推断闭合和开放状态，DNA发夹两端连接有荧光供体F和猝灭剂Q

新标签中打开下载幻灯片

在活细胞中，随着两条DNA链（双螺旋）之间分子间配对的断裂，松散的链通常会形成DNA发夹结构。DNA发夹结构参与许多重要的生物学功能，包括基因表达的调节（Zazopoulos等。,1997)、DNA重组（Froelich-Ammon等。,1994)以及诱变事件的促进作用（Trinh和Sinden，1993). 发夹结构也可能是潜在的反义药物（Tang等。,1993)因为将发夹注射到活细胞中，发夹的核酸碱基与疾病基因的信使核糖核酸互补，可能会阻止该基因的表达。由于其生物学相关性，研究DNA发夹的构象特性，如开态和闭态之间的构象涨落和能量屏障，是理解更复杂生化过程的重要模型系统。

1.2. 荧光寿命实验

DNA发夹自发地在开放和闭合状态之间切换。在正在研究的单分子荧光寿命实验中，一个荧光供体染料和一个猝灭剂附着在DNA分子的两端(图1). 然后，将分子放置在激光脉冲列照亮的焦点体积中。施主染料被激光脉冲激发时会发射光子，而猝灭剂会消除激发。在发夹的闭合状态下，猝灭很强，很少检测到来自施主的光子；在开放态下，由于施主和猝灭体相距较远，猝灭较弱，许多来自施主的光子被检测到。对于检测到的光子，将记录其到达时间。通过将检测到的光子到达时间与常规激光脉冲列进行比较，我们可以进一步确定施主染料从激发时刻释放光子所需的时间长度。这个时间长度称为光子延迟时间。图2说明了到达时间的确定t吨和延迟时间τ。

激光激发供体染料，从而释放光子：通过将检测到的光子到达时间与常规激光脉冲串进行比较，确定光子延迟时间（相对于其激发脉冲）

图2

激光激发施主染料，施主染料释放光子：通过比较检测到的光子到达时间和常规激光脉冲序列，确定光子延迟时间（相对于其激发脉冲）

新标签中打开下载幻灯片

让A类表示闭合，并让B类表示打开。DNA发夹动力学最简单的模型是一个连续时间的两态马尔可夫链，称为二态模型（雷利和斯金纳，1993,1994年a,b条)：从状态开始的分子A类会留在家里A类等待时间呈指数分布，然后跳到状态B类; 反之亦然B类。这可以描述为

A类 \underset{{k个}_{21}}{\overset{{k个}_{12}}{⇄}} B类,

(1.1)

哪里k个₁₂和k个₂₁表示两个状态之间的转换速率常数。让

问 = (\begin{matrix} - {k个}_{12} & {k个}_{12} \\ {k个}_{21} & - {k个}_{21} \end{matrix})

是强度矩阵（无穷小生成器）。Kolmogorov正向方程（Fokker–Planck方程）给出了过渡矩阵

P（P） (t吨) = 经验 (问 t吨) = (\begin{matrix} π_{1} + π_{2} 经验 (- k个 t吨) & π_{2} {1 - 经验 (- k个 t吨)} \\ π_{1} {1 - 经验 (- k个 t吨)} & π_{2} + π_{1} 经验 (- k个 t吨) \end{matrix}) 。

(1.2)

哪里k个=k个₁₂+k个₂₁和(π₁,π₂)=(k个₂₁/(k个₁₂+k个₂₁),k个₁₂/(k个₁₂+k个₂₁))是两状态马尔可夫链的稳态分布。

光子到达时间t吨以及相应的延迟时间τ取决于底层隐藏的双状态进程。让γ(t吨)表示DNA发夹结构的荧光水平t吨，它接受值一和b条分别以州为单位A类和B类具有一>b条光子到达时间t吨遵循双随机泊松过程，到达率与γ(t吨)其中，由于到达时间和到达率都是随机的，因此出现了术语“双重随机”。在文献中，此过程也称为Cox过程（Cox，1955). 对于检测到的光子，其光子延迟时间τ遵循指数分布γ(t吨). 延迟时间对基础状态的依赖性是因为，当施主和猝灭体彼此靠近时（即在闭合状态下），能量从施主转移到猝灭体，从而缩短延迟时间。图3说明了数据对两个状态进程的依赖性γ粗体线描述了与DNA发夹的开-闭状态相对应的未观察到的二态马尔可夫链。每个垂直钢筋图3表示光子的到达（时间）。我们可以看到，在打开状态下（即γ⁻¹)光子到达的频率更高。对于每个条形图（即每个光子到达），高度表示相应的延迟时间τ记录的光子。在打开状态下，延迟时间τ往往更大。因此，DNA发夹的打开或关闭状态可以从检测到的光子到达时间和延迟时间（Lu等。,1998; 贾等。,1997; 埃格林等。,1998).

单分子寿命实验中的数据结构：|，光子到达时间t，单位为秒（高度表示延迟时间τ，单位为纳秒），未观测到的两状态马尔可夫链，t和τ都依赖于此

图3

单分子寿命实验中的数据结构：|，光子到达时间t吨以秒为单位（高度表示延迟时间τ以纳秒为单位）；，未观察到的两状态马尔可夫链，其中t吨和τ依赖

新标签中打开下载幻灯片

除了隐马尔可夫结构之外未观察到的分子在溶剂中的布朗扩散轨迹增加了另一层复杂性，使推断大大复杂化（参见第节三). 为了有效地解决这些问题，我们使用贝叶斯数据增强方法，“插补”缺失的布朗过程（Tanner和Wong，1987).

随着越来越多的分子数据的积累，科学家开始质疑简单的两态模型。例如，有人认为，两个状态之间的能量屏障也可以动态波动，或者两个状态中的每一个都可能有子状态，并且这些子状态可以以不同的速率进行通信（Cao，2000; 杨和曹，2001,2002). 根据目前的数据分辨率和现有的推断方法，识别不同的模型并评估其与实验数据的拟合程度仍然难以捉摸，无法仅仅猜测（见Schenter等。(1999)以供进一步讨论）。使用贝叶斯数据增强方法，我们能够提供重要的统计证据，从而得出以下结论：简单的两状态模型不允许能垒在两种状态之间波动，不足以解释数据。揭示能量屏障的本质可能是全面了解DNA发夹生物物理特性的第一步。

本文组织如下。章节2详细介绍了基于闭合形式似然的两状态统计模型和相应的贝叶斯分析。章节三介绍了处理实验复杂性（如分子布朗运动）的数据增强方法。章节4考虑两个国家以外的模式。章节5分析实验和模拟数据并讨论模型评估问题，其中，通过研究连续扩散模型，我们证明简单的两态模型不足以捕捉DNA发夹构象动力学的精细细节。章节6以进一步的讨论结束本文。

2.两国模型的贝叶斯分析

让Y（Y）(t吨)是最新到达的光子总数t吨.在无限小的时间间隔内(t吨,t吨+d日t吨)观察光子的概率与γ⁻¹(t吨)d日t吨.表示Y（Y）_t吨=Y（Y）(t吨+d日t吨)−Y（Y）(t吨)，我们有

P（P） (Δ {Y（Y）}_{t吨} = 1 | γ_{t吨}) = {A类}_{0} (t吨) γ^{- 1} (t吨) d日 t吨,

(2.1)

[τ | Δ {Y（Y）}_{t吨} = 1, γ_{t吨}] ~ γ (t吨) 经验 {- γ (t吨) τ}

(2.2)

哪里A类₀(t吨)是光子到达时的强度t吨，可能会随时间变化。随机性A类₀(t吨)将在第节中说明三在本节中，为了集中讨论基本思想，我们将A类₀(t吨)作为随时间变化的常数A类₀(t吨)≡A类₀按照随机过程文献中的约定，我们使用括号和下标来索引（相同的）随机过程；例如γ(t吨)和γ_t吨两者都表示底层的两状态马尔可夫链，以及Y（Y）(t吨)和Y（Y）_t吨两者都涉及光子计数过程。重要的是要指出这个公式(2.1)等效于Cox过程中的光子到达，但使用方程式（2.1）我们将演示一个简单的概括。

2.1. 可能性计算

让0=t吨₀<t吨₁< …<t吨_n个是观测到的光子到达时间，让τ_我是相应的光子延迟时间。这对搭档 ${({t吨}_{我}, τ_{我})}_{我 = 0}^{n个}$ 通过荧光寿命实验收集。请注意，可能性由三部分组成：

（a）
光子到达时间t吨_我;
（b）
延迟时间τ_我;
（c）
没有光子以时间间隔到达的事实(t吨_我,t吨_我+1).

还要注意，Cox过程可以看作是一个隐马尔可夫过程。但它不同于标准的隐马尔可夫模型（Elliott等。,1999)相邻光子到达之间的间隔长度也是隐藏过程的信息γ(t吨). 鉴于γ=(γ(t吨₀),…,γ(t吨_n个))，联合似然函数可以写成

\begin{matrix} L（左） (t吨, τ, γ ∣ θ) = P（P） {γ ({t吨}_{0})} [\prod_{我 = 0}^{n个} {\frac{{A类}_{0}}{γ ({t吨}_{我})} d日 t吨}] [\prod_{我 = 0}^{n个} γ ({t吨}_{我}) 经验 {- γ ({t吨}_{我}) τ_{我}}] \\ \times \prod_{我 = 0}^{n个 - 1} P（P） {Y（Y） ({t吨}_{我 + 1}^{-}) - Y（Y） ({t吨}_{我}) = 0, γ ({t吨}_{我 + 1}) ∣ γ ({t吨}_{我})} 。 \end{matrix}

使用无穷小矩阵方法，在附录A我们证明了这一点

P（P） {Y（Y） ({t吨}_{我 + 1}^{-}) - Y（Y） ({t吨}_{我}) = 0, γ ({t吨}_{我 + 1}) ∣ γ ({t吨}_{我})} = {[经验 {(问 - H（H）) ({t吨}_{我 + 1} - {t吨}_{我})}]}_{(我 ({t吨}_{我}), 我 ({t吨}_{我 + 1})),}

哪里问如第节所示1.2,

H（H） = (\begin{matrix} {A类}_{0} / 一 & 0 \\ 0 & {A类}_{0} / b条 \end{matrix})

和索引我(t吨)等于1，如果γ(t吨) =一否则为2。因此，无条件似然可以通过假设γ(t吨₀)从平稳分布开始，求出γ(t吨_我)通过递归矩阵乘法：

L（左） (t吨, τ ∣ θ) = (π_{1}, π_{2}) {D类}_{0} H（H） [\prod_{我 = 0}^{n个 - 1} 经验 {(问 - H（H）) Δ {t吨}_{我}} {D类}_{我 + 1} H（H）] (\begin{array}{l} 1 \\ 1 \end{array})

(2.3)

其中(π₁,π₂)是的平衡分布γ(t吨), Δt吨_我=t吨_我+1−t吨_我和

{D类}_{我} = (\begin{matrix} 一 经验 (- 一 τ_{我}) & 0 \\ 0 & b条 经验 (- b条 τ_{我}) \end{matrix}) 。

备注1。我们注意到，通过第一步分析，也可以从无穷小生成器开始，建立一组微分方程（Karlin和Taylor，1998; 卡尔，1986). 与第一步分析相比附录A更直观，更容易概括，以适应其他复杂情况，如二乘二模型（Schenter等。,1999)和分子布朗扩散（参见第节三). 例如，在二乘二模型中，假设两个状态中的每一个都有两个子状态：

\begin{array}{l} {A类}_{1} & ⇄_{{k个}_{21}}^{{k个}_{12}} & {B类}_{1} \\ α ↓ ↑ β & α ↓ ↑ β \\ {A类}_{2} & ⇄_{{k个}_{21}^{'}}^{{k个}_{12}^{'}} & {B类}_{2} \end{array}

(2.4)

哪里α和β是子状态之间的转换速率。这个γ(t吨)过程有价值一₁,一₂,b条₁和b条₂在美国A类₁,A类₂,B类₁和B类₂分别是。相同的似然公式(2.3)与保持问,H（H）和D类_我被新模型中的相应矩阵替换：

\begin{matrix} 问 = (\begin{array}{l} - ({k个}_{12} + α) & α & {k个}_{12} & 0 \\ β & - ({k个}_{12}^{'} + β) & 0 & {k个}_{12}^{'} \\ {k个}_{21} & 0 & - ({k个}_{21} + α) & α \\ 0 & {k个}_{21}^{'} & β & - ({k个}_{21}^{'} + β) \end{array}) \begin{array}{l} \leftarrow {A类}_{1} \\ \leftarrow {A类}_{2} \\ \leftarrow {B类}_{1} \\ \leftarrow {B类}_{2} \end{array} \\ H（H） = 诊断 (\frac{{A类}_{0}}{一_{1}}, \frac{{A类}_{0}}{一_{2}}, \frac{{A类}_{0}}{{b条}_{1}}, \frac{{A类}_{0}}{{b条}_{2}}), \\ {D类}_{我} = 诊断 {一_{1} 经验 (- 一_{1} τ_{我}), 一_{2} 经验 (- 一_{2} τ_{我}), {b条}_{1} 经验 (- {b条}_{1} τ_{我}), {b条}_{2} 经验 (- {b条}_{2} τ_{我})} 。 \end{matrix}

我们将关于二乘二模型的进一步讨论推迟到第节4，其中考虑了两个国家模型以外的模型。值得注意的是，Wolpert和Ickstadt(1998)已将Cox过程推广到双随机层次模型，以考虑点过程模型中基本速率的不确定性和空间变化。

2.2. 马尔可夫链蒙特卡罗算法

两国模式有五个参数：一,b条,π₁,k个和A类₀.让η(θ)表示参数的先验分布θ=(一,b条,π₁,k个,A类₀). 那么后验分布是

P（P） (θ ∣ t吨, τ) \propto η (θ) L（左） (t吨, τ ∣ θ),

其中包含以下约束：

（a）
一>b条> 0,
（b）
0 ≼π₁≼1，
（c）
k个> 0,
（d）
A类₀> 0.

由于直接计算自P（P）(θ|t吨,τ)是不可行的，我们设计了一个Metropolis-type算法来从中进行采样。该算法迭代了以下步骤。

第1步：给定一和b条,

（i）
画x个来自Γ（1/c（c）₁,交流电₁)和年来自Γ（1/c（c）₂,公元前₂)，其中c（c）₁和c（c）₂是两个调谐参数，以及
（ii）
让一^′=最大值(x个,年)和b条^′=最小值(x个,年).

调谐参数c（c）₁和c（c）₂控制提案的步长。拿一^′和b条^′分别为最大值和最小值x个和年满足约束（a），从而使建议密度成为两个产品伽马密度的混合物。

第2步：给定π₁，绘制 $π_{1}^{'}$ 从beta分布B类{c（c）_三π₁,c（c）_三(1−π₁)}. 由于B类{c（c）_三π₁,c（c）_三(1−π₁)}是π₁和π₁(1−π₁)/(c（c）_三+1）分别是， $π_{1}^{'}$ 可视为局部扰动π₁、和调整参数c（c）_三控制扰动的局部性。

步骤3：给定k个，绘制k个^′来自Γ（1/c（c）₄,千卡₄). 平均值和方差k个^′/k个是1和c（c）₄分别出租c（c）₄微调扰动。

第4步：给定A类₀，绘制 ${A类}_{0}^{'}$ 来自Γ（1/c（c）₅,A类₀c（c）₅)，其平均值为A类₀，其方差由c（c）₅

该算法在我们的815个数据对的模拟研究中运行良好 ${({t吨}_{我}, τ_{我})}_{我 = 0}^{815}$ 从隐藏的γ过程（815是实际实验中典型的数据大小）。图4总结了参数的后验抽样分布（具有平坦的先验），其中竖线代表真实值。该算法除了正确识别所有参数外，运行速度也很快：在奔腾4个人计算机上绘制10000个样本所用的时间不到2分钟图4.快速计算主要是因为直接计算似然函数(2.3)只需要O（运行）(n个)操作。

模拟二态数据的后验样本直方图（，真值）：（a）a；（b） b；（c） π1；（d） k；（e） A0级

图4

模拟两国数据的后验样本直方图(，真值）：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀

新标签中打开下载幻灯片

3.两态模型中的布朗扩散

3.1. 统计公式

在单分子寿命实验中，当激光激发的染料发射光子时，DNA发夹分子也会在焦点体积中扩散。因此，由于DNA分子接收到的激光照明强度在焦体积中心最高，并从中心向外减小，因此实际光子到达强度A类₀(t吨)变化。数学上，我们可以写A类₀(t吨) =A类₀ α(t吨)带有

α (t吨) = 经验 {- \frac{{B类}_{x个}^{2} (t吨) + {B类}_{年}^{2} (t吨)}{2 {w个}_{x个 年}^{2}} - \frac{{B类}_{z（z）}^{2} (t吨)}{2 {w个}_{z（z）}^{2}}},

(3.1)

其中(B类_x个(t吨),B类_年(t吨)和B类_z（z）(t吨))是DNA分子的位置和常数w个_xy公司和w个_z（z）指定x个-,年-和z（z）-椭球焦点体积的轴。让σ表示扩散常数，即。B类_x个(t吨) =σ W公司₁(t吨),B类_年(t吨) =σ W公司₂(t吨)和B类_z（z）(t吨) =σ W公司_三(t吨)，其中(W公司₁(t吨),W公司₂(t吨),W公司_三(t吨))是三个独立的标准布朗运动。在寿命实验中，常数w个_xy公司,w个_z（z）和σ根据荧光相关光谱理论（Magde等。,1974)并假定已知。在存在扩散的情况下，概率公式(2.1)和(2.2)更改为

P（P） (Δ {Y（Y）}_{t吨} = 1 ∣ γ_{t吨}, α_{t吨}) = {A类}_{0} α (t吨) γ^{- 1} (t吨) d日 t吨,

（3.2）

[τ ∣ Δ {Y（Y）}_{t吨} = 1, γ_{t吨}, α_{t吨}] \sim γ (t吨) 经验 {- γ (t吨) τ} 。

(3.3)

由于对强度的额外限制α(t吨)，可能性也会发生变化。因为在实验中α(t吨)是相对于光子到达的一个慢变过程，我们做了一个近似α(t吨)≈α(t吨_我)的t吨∈ (t吨_我,t吨_我+1). 然后，使用与中相同的离散化技术附录A，我们有

P（P） {Y（Y） ({t吨}_{1}^{-}) - Y（Y） ({t吨}_{0}) = 0, γ ({t吨}_{1}) ∣ γ ({t吨}_{0}), α ({t吨}_{0})} = {(经验 [{问 - α ({t吨}_{0}) H（H）} ({t吨}_{1} - {t吨}_{0})])}_{(我_{{t吨}_{0}}, 我_{{t吨}_{1}})},

其中矩阵问和H（H）和以前一样。类似于节2，结合区间的似然贡献{(t吨_我,t吨_我+1);我= 0,1,…,n个−1}与配对{t吨_我,τ_我}产生给定的条件似然α_t吨，修改了方程式（2.3）以下为：

L（左） (t吨, τ ∣ θ, α_{t吨}) = (π_{1}, π_{2}) {D类}_{0} {H（H）}_{0} [π_{我 = 0}^{n个 - 1} 经验 {(问 - {H（H）}_{我}) Δ {t吨}_{我}} {D类}_{我 + 1} {H（H）}_{我 + 1}] (\begin{array}{l} 1 \\ 1 \end{array}),

(3.4)

哪里H（H）_我表示α(t吨_我)H（H）=α(t吨_我)诊断(A类₀/一,A类₀/b条). 事先分配η(θ)关于参数的后验分布θ根据观察结果 ${({t吨}_{我}, τ_{我})}_{我 = 0}^{n个}$ 是

P（P） (θ ∣ t吨, τ) \propto \int η (θ) L（左） (t吨, τ ∣ θ, α_{t吨}) P（P） (α_{t吨}) d日 (α_{t吨}),

(3.5)

哪里P（P）(一_t吨)表示的概率定律α(t吨). 比较表达式(3.5)节中有对应的表达式2，我们注意到布朗扩散因子α(t吨)分子的不可观察性，必须整合出来。该路径积分无法进行分析计算，必须通过数据增强方法进行处理（Tanner和Wong，1987).

3.2. 数据增强

直观地说，如果布朗扩散(B类_x个,B类_年,B类_z（z）)给出了第节中描述的计算2可以很容易地执行。因此，通过从关节后部分布抽样θ和(B类_x个,B类_年,B类_z（z）),

P（P） (θ, {B类}_{x个}, {B类}_{年}, {B类}_{z（z）} ∣ t吨, τ) \propto η (θ) L（左） (t吨, τ ∣ θ, α_{t吨}) P（P） ({B类}_{x个}) P（P） ({B类}_{年}) P（P） ({B类}_{z（z）}),

（3.6）

其中的依赖性α_t吨上的(B类_x个,B类_年,B类_z（z）)由提供方程式（3.1），我们有效地边缘化了隐藏的扩散过程，并获得了正确的推断。我们只需要在光子到达时间和相应的每个表达式项处增加三维布朗运动(3.6)具有简单形式。例如B类_x个是

P（P） ({B类}_{x个}) = {(2 π)}^{- n个 / 2} 经验 [- \frac{1}{2 σ^{2}} \sum_{我 = 0}^{n个 - 1} \frac{{{B类}_{x个} ({t吨}_{我 + 1}) - {B类}_{x个} ({t吨}_{我})}^{2}}{Δ {t吨}_{我}}] 。

(3.7)

从首字母开始θ和配置(B类_x个,B类_年,B类_z（z）)，我们迭代以下条件采样步骤以从表达式模拟(3.6)。

第1步：绘制θ电流扩散条件(B类_x个,B类_年,B类_z（z）),

[θ ∣ {B类}_{x个}, {B类}_{年}, {B类}_{z（z）}, t吨, τ] \sim P（P） (θ ∣ {B类}_{x个}, {B类}_{年}, {B类}_{z（z）}, t吨, τ) \propto η (θ) L（左） (t吨, τ ∣ θ, α_{t吨}) 。

(3.8)

第2步：绘制扩散(B类_x个,B类_年,B类_z（z）)调节电流值θ,

[{B类}_{x个}, {B类}_{年}, {B类}_{z（z）} ∣ θ, t吨, τ] \sim P（P） ({B类}_{x个}, {B类}_{年}, {B类}_{z（z）} ∣ θ, t吨, τ) \propto L（左） (t吨, τ ∣ θ, α_{t吨}) P（P） ({B类}_{x个}) P（P） ({B类}_{年}) P（P） ({B类}_{z（z）}) 。

(3.9)

第一步由Metropolis–Hastings（MH）算法实现，如第节所述2由于MH算法中涉及的计算比步骤2中的采样简单得多，特别是当链较长时，我们在步骤1中多次进行MH迭代，然后再转到步骤2。

步骤2是通过逐个时间点更新扩散链时间点来实现的。换句话说，对于每个我= 0,1,…,n个，我们建议一个新位置 $\tilde{B类} = (\tilde{x个}, \tilde{年}, \tilde{z（z）})$ 对于我第个时间点B类(t吨_我)=(B类_x个(t吨_我),B类_年(t吨_我),B类_z（z）(t吨_我))并根据MH规则接受建议。由于每次只考虑一个时间点，因此MH比率很容易计算。例如，如果我第个时间点B类(t吨_我)=(B类_x个(t吨_我),B类_年(t吨_我),B类_z（z）(t吨_我))建议更改为 $\tilde{B类} = (\tilde{x个}, \tilde{年}, \tilde{z（z）})$ ⁠，MH概率很简单

最小值 [1, \frac{L（左） (t吨, τ ∣ θ, {\tilde{α}}_{t吨})}{L（左） (t吨, τ ∣ θ, α_{t吨})} \frac{P（P） ({\tilde{B类}}_{x个}) P（P） ({\tilde{B类}}_{年}) P（P） ({\tilde{B类}}_{z（z）})}{P（P） ({B类}_{x个}) P（P） ({B类}_{年}) P（P） ({B类}_{z（z）})} \frac{T型 {\tilde{B类} \to B类 ({t吨}_{我})}}{T型 {B类 ({t吨}_{我}) \to \tilde{B类}}}],

哪里T型{·→·}是提案的过渡密度

\frac{P（P） ({\tilde{B类}}_{x个}) P（P） ({\tilde{B类}}_{年}) P（P） ({\tilde{B类}}_{z（z）})}{P（P） ({B类}_{x个}) P（P） ({B类}_{年}) P（P） ({B类}_{z（z）})}

以下为方程式（3.7）.评估似然函数L（左）(t吨,τ∣θ,α_t吨)我们高效地实现了以下内容向前–向后递归方案：反向计算似然中的部分和（积分）；一次向前更新一个组件（时间点）的扩散链。这种方法确保计算成本仅为O（运行）(n个)操作，而天真的方法是有序的O（运行）(n个²).

为了开始前向-后向算法，我们首先计算向后的矩阵序列K（K）_我通过递归

{\begin{array}{l} {K（K）}_{n个 + 1} = 我, {K（K）}_{n个} = {D类}_{n个} {H（H）}_{n个}, \\ {K（K）}_{我} = {D类}_{我} {H（H）}_{我} 经验 {(问 - {H（H）}_{我}) Δ {t吨}_{我}} {K（K）}_{我 + 1}, 我 < n个, \end{array}

其中，与之前一样，D类_我=诊断{一经验（−一τ_我),b条经验（−b条τ_我)},H（H）_我=α(t吨_我)诊断(A类₀/一,A类₀/b条)和Δt吨_我=t吨_我+1−t吨_我.表示v（v）₀=(π₁,π₂). 然后向前地采样算法可以描述如下我= 0,1,…,n个。

第1步：提出更改 $\tilde{B类} = (\tilde{x个}, \tilde{年}, \tilde{z（z）})$ 对于我第个时间点B类(t吨_我).

第2步：计算

对 = {\begin{array}{l} {D类}_{我} {H（H）}_{我} 经验 {(问 - {H（H）}_{我}) Δ {t吨}_{我}} & 如果 我 < n个, \\ {D类}_{n个} {H（H）}_{n个} & 如果 我 = n个, \end{array}

和

S公司 = {\begin{array}{l} {D类}_{我} {\tilde{H（H）}}_{我} 经验 {(问 - {\tilde{H（H）}}_{我}) Δ {t吨}_{我}} & 如果 我 < n个, \\ {D类}_{n个} {\tilde{H（H）}}_{n个} & 如果 我 = n个, \end{array}

哪里 ${\tilde{H（H）}}_{我} = \tilde{α} ({t吨}_{我}) d日我一克 ({A类}_{0} / 一, {A类}_{0} / b条)$ 和

\tilde{α} ({t吨}_{我}) = 经验 (- \frac{{\tilde{x个}}^{2} + {\tilde{年}}^{2}}{2 {w个}_{x个 年}^{2}} - \frac{{\tilde{z（z）}}^{2}}{2 {w个}_{z（z）}^{2}})

是照明强度的值我提案的第个时间点 $\tilde{B类}$ ⁠因此， ${v（v）}_{我} 对 {K（K）}_{我 + 1} (\begin{array}{l} 1 \\ 1 \end{array})$ 现在是可能的L（左）(t吨,τ|θ,α_t吨)在原始扩散链上进行评估，以及 ${v（v）}_{我} S公司 {K（K）}_{我 + 1} (\begin{array}{l} 1 \\ 1 \end{array})$ 是可能性 $L（左） (t吨, τ ∣ θ, {\tilde{α}}_{t吨})$ 在扩散链上用新提出的我第个时间点。向量v（v）_我用于记录第一个我时间点；它将在步骤4中被更新。

步骤3：计算MH比率

第页 = \frac{L（左） (t吨, τ ∣ θ, {\tilde{α}}_{t吨}) P（P） ({\tilde{B类}}_{x个}) P（P） ({\tilde{B类}}_{年}) P（P） ({\tilde{B类}}_{z（z）}) T型 {\tilde{B类} \to B类 ({t吨}_{我})}}{L（左） (t吨, τ ∣ θ, α_{t吨}) P（P） ({B类}_{x个}) P（P） ({B类}_{年}) P（P） ({B类}_{z（z）}) T型 {B类 ({t吨}_{我}) \to \tilde{B类}}} 。

第4步：生成u个∼U型(0,1). 如果u个<最小值（1，第页)，我们更新B类(t吨_我)至 $\tilde{B类}$ 然后让v（v）_我+1=v（v）_我S公司; 否则，我们会继续B类(t吨_我)不变并出租v（v）_我+1=v（v）_我对.矢量v（v）_我+1现在记录第一个我+1个时间点。

注意，由于直接（有条件）采样B类(t吨_我)很难。虽然可以对B类(t吨)，即更新流程的一部分(B类(t吨_我),B类(t吨_我+1),…,B类(t吨_我+b条−1))同时，由于多维MH算法难以实现，因此不清楚这是否会显著改善收敛性。此外，我们发现上述分量更新方案利用了有效的似然评估方案，并且已经非常有效。在后面的部分中，我们将说明分块移动对于处理更复杂的模型非常有用。

3.3. 其他实验并发症

由于技术限制，除了分子布朗扩散外，实际实验还有更多的复杂性。第一个额外的复杂性来自背景光子。到目前为止，在分析中，我们假设所记录的光子都来自激光激发的DNA发夹结构。但在实际实验中，激光也会从背景中激发光子，光子的到达被描述为具有速率的时间均匀泊松过程ρ.参数ρ可以通过运行一个长时间的控制实验来估计，在焦距中没有分子（因此检测到的光子都来自背景），并且我们假定已知。因此，我们只需要从以下方面改变整个光子到达过程（Cox过程）方程式（3.2）到

P（P） (Δ {Y（Y）}_{t吨} = 1 ∣ γ_{t吨}, α_{t吨}) = {ρ + \frac{{A类}_{0}}{γ (t吨)} α_{t吨}} d日 t吨 。

(3.10)

似然函数从方程式（3.4）到

L（左） (t吨, τ ∣ θ, α_{t吨}) = (π_{1}, π_{2}) {D类}_{0} {G公司}_{0} [\prod_{我 = 0}^{n个 - 1} 经验 {(问 - {G公司}_{我}) Δ {t吨}_{我}} {D类}_{我 + 1} {G公司}_{我 + 1}] (\begin{array}{l} 1 \\ 1 \end{array}),

（3.11）

其中矩阵D类_我=诊断{一经验（−一τ_我),b条经验（−b条τ_我)}如前所述，并且新矩阵

{G公司}_{我} = d日 我 一 克 {ρ + \frac{{A类}_{0}}{一} α ({t吨}_{我}), ρ + \frac{{A类}_{0}}{b条} α ({t吨}_{我})} 。

(3.12)

第二个实验复杂性是与光子延迟时间测定相关的时间缠绕问题。如第节所示1.2延迟时间是通过比较检测到的光子到达时间和常规激光脉冲序列来确定的。然而，仔细检查图2这表明，如果延迟时间超过两个相邻激光脉冲之间的间隔13.2ns，我们就不再确定是哪个激光脉冲激发了供体染料。机器读数假设染料被检测到的光子到达之前的脉冲激发，但这会产生时间缠绕效应。例如，如果实际延迟时间为17 ns，机器将确定为17−13.2=3.8 ns。第三个复杂因素来自负读数：一些光子延迟时间τ记录为负值。负面消息τ表示光子已到达，但缺少延迟时间。通过进一步修改似然函数，可以相对容易地处理这两种复杂情况。由于指数分布的无记忆特性，时间包装改变了延迟时间的分布τ转换为截断的指数分布，从而改变表达式(3.3)到

[τ ∣ Δ {Y（Y）}_{t吨} = 1, γ_{t吨}, α_{t吨}] \sim \frac{γ (t吨) 经验 {- γ (t吨) τ}}{1 - 经验 {- γ (t吨) M（M）}},

(3.13)

哪里M（M）是包裹点（上例中为13.2 ns）。负片记录进一步改变了分布(3.13)到

[τ ∣ Δ {Y（Y）}_{t吨} = 1, γ_{t吨}, α_{t吨}, τ > 0] \sim \frac{γ (t吨) 经验 {- γ (t吨) τ}}{1 - 经验 {- γ (t吨) M（M）}} 。

因此，更换D类_我-矩阵依据

{\tilde{D类}}_{我} = {\begin{array}{l} d日 我 一 克 {\frac{一 经验 (- 一 τ_{我})}{1 - 经验 (- 一 M（M）)}, \frac{b条 经验 (- b条 τ_{我})}{1 - 经验 (- b条 M（M）)}} & 如果 τ_{我} ⩾ 0, \\ 我 & 如果 τ_{我} < 0 \end{array}

(3.14)

在公式中（3.11）给出了真实实验数据的可能性。

4.超越两态模型：连续扩散模型

4.1. 模型

两国模式(1.1)将构象动力学描述为两状态连续时间马尔可夫链。然而，科学家观察到，对于某些蛋白质和酶，两态模型不足以准确描述构象细节。例如，两态模型暗示了不同状态下连续等待时间之间的“无记忆”特性，而在许多涉及蛋白质和酶的单分子实验中，这种特性被认为是违反的（Schenter等。,1999; 卢等。,1998). 这种现象被称为“动力障碍”（谢，2002; 茨万齐格，1990)，推动了超越两国模式的模式。两个接两个模型（2.4）第节简要介绍了2就是这样的尝试。在这个模型中，状态A类₁和A类₂实验上无法区分B类₁和B类₂。因此，从实验中可以观察到的是A类₁和A类₂（和分别B类₁和B类₂). 在统计学文献中，这类过程被称为聚集马尔可夫过程（Fredkin和Rice，1986). 聚集的结果是一个分子在集体态中所花费的时间{A类₁,A类₂}和{B类₁,B类₂}不再独立，这可以解释动力紊乱（曹，2000; 杨和曹，2001,2002). 如第节所示2，二乘二模型的似然函数与两国模型的形式相同。因此，第节中概述的贝叶斯数据增强方法2和三适用于两乘二模型。推理过程基本上是相同的，只是二乘二模型中有更多的参数。

二乘二模型可以进一步推广为连续扩散模型，该模型无法通过早期开发的数据增强方法直接处理，这也是本节的重点。在这个模型中，不是使用离散聚合来解释动态无序，而是一个连续的随机控制过程x个(t吨)引入了“控制”转换率的方法，如下所示：

A类 \underset{{k个}_{21} 经验 {- x个 (t吨)}}{^{{k个}_{12} 经验 {- x个 (t吨)}}} \underset{}{⇄} B类 。

(4.1)

从物理上看，这可以看作是两种状态之间能量屏障动态波动的结果(x个(t吨)与能量屏障有关小时(t吨)通过x个(t吨) =小时(t吨)/k个_B类T型，其中k个_B类是玻尔兹曼常数T型是温度）。这个波动的能量屏障（参见图5举例来说）允许转换速率是时变的和随机的。相应地，各州的连续等待时间A类和B类不是独立的。

两种状态之间的能垒波动：奥恩斯坦-乌伦贝克过程模拟了DNA发夹打开和关闭状态之间能垒的动态振荡

图5

两种状态之间的涨落势垒：Ornstein–Uhlenbeck过程模拟了DNA发夹打开和关闭状态之间能量势垒的动态振荡

新标签中打开下载幻灯片

正如Agmon和Hopfield首先建议的那样(1983)，Ornstein–Uhlenbeck过程用于捕捉两个状态之间能量屏障的波动：

d日 {x个}_{t吨} = - λ {x个}_{t吨} d日 t吨 + \sqrt{(2 ξ λ) d日 {W公司}_{t吨} 。}

(4.2)

一种描绘的方式模型（4.1）就是认为DNA分子具有无限多的无法区分的状态{A类₁,A类₂、…}和{B类₁,B类₂，…}，二乘二模型的连续推广。尽管已经开始了许多辩论（博内等。,1998; 应等。,2001; 格伦威尔等。,2001; 安萨里等。,2001,2002)到目前为止，对于DNA发夹的连续扩散模型是否确实比简单的两态模型更合适，还没有明确的证据。因此，最有效地综合所有可用的统计证据以区分模型是很有意义的。

4.2. 连续扩散模型下的似然性和数据扩充

采用控制过程的观察结果x个(t吨)在两个连续光子到达之间是稳定的，我们采用近似值x个(t吨)≈x个(t吨_我)的t吨∈ (t吨_我,t吨_我+1)并应用第节中的技术2获得闭合形式的条件似然α(t吨)和x个(t吨):

L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) = (π_{1}, π_{2}) {\tilde{D类}}_{0} {G公司}_{0} [\prod_{我 = 0}^{n个 - 1} 经验 {(问_{我} - {G公司}_{我}) Δ {t吨}_{我}} {\tilde{D类}}_{我 + 1} {G公司}_{我 + 1}] (\begin{array}{l} 1 \\ 1 \end{array}),

(4.3)

哪里G公司_我和 ${\bar{D类}}_{我}$ 如前所述方程式（3.12）和(3.14)、和

问_{我} ≜ (\begin{matrix} - {k个}_{12} 经验 {- x个 ({t吨}_{我})} & {k个}_{12} 经验 {- x个 ({t吨}_{我})} \\ {k个}_{21} 经验 {- x个 ({t吨}_{我})} & - {k个}_{21} 经验 {- x个 ({t吨}_{我})} \end{matrix}) 。

参数的后验分布(θ,λ,ξ)根据观察结果{(t吨_我,τ_我)}是

P（P） (θ, λ, ξ ∣ t吨, τ) \propto \iint η (θ, λ, ξ) L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） (α_{t吨}) P（P） ({x个}_{t吨} ∣ λ, ξ) d日 (α_{t吨}) d日 ({x个}_{t吨}),

(4.4)

哪里P（P）(x个_t吨|λ,ξ)是过程的过渡密度x个(t吨)和η表示参数的先验分布。

通过增强两者x个(t吨)和布朗扩散B类(t吨)，我们考虑(θ,λ,ξ)和(B类(t吨),x个(t吨)):

P（P） (θ, λ, ξ, B类, x个 ∣ t吨, τ) \propto η (θ, λ, ξ) L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） {B类 (t吨)} P（P） ({x个}_{t吨} ∣ λ, ξ) 。

(4.5)

从Kolmogorov反向方程

\frac{\partial P（P）}{\partial t吨} = ξ λ \frac{\partial^{2} P（P）}{\partial {x个}^{2}} - λ x个 \frac{\partial P（P）}{\partial x个},

我们有x个(t吨_我+1)|x个(t吨_我)∼N个[x个(t吨_我)pt（磅）经验（−λΔt吨_我),ξ{1−pt（磅）经验（-2λΔt吨_我)}]（见Karlin和Taylor(1981)，第218页）。假设Ornstein–Uhlenbeck过程x个(t吨)从平稳分布开始，我们得到了联合密度

\begin{array}{l} P（P） ({x个}_{t吨} ∣ λ, ξ) \propto ξ^{- (n个 + 1) / 2} [\prod_{我 = 0}^{n个 - 1} {1 - 经验 (- 2 λ Δ {t吨}_{我})}^{- 1 / 2}] \\ \times 经验 [- \frac{x个 {({t吨}_{0})}^{2}}{2 ξ} - \sum_{我 = 0}^{n个 - 1} \frac{{x个 ({t吨}_{我 + 1}) - x个 ({t吨}_{我}) 经验 (- λ Δ {t吨}_{我})}^{2}}{2 ξ {1 - 经验 (- 2 λ Δ {t吨}_{我})}}] 。 \end{array}

(4.6)

为了使马尔可夫链蒙特卡罗抽样更有效，我们重新参数化了控制过程x个(t吨)通过更换ξ具有

ϕ = \sqrt{(ξ λ)} 。

然后将Ornstein–Uhlenbeck过程更改为

d日 {x个}_{t吨} = - λ {x个}_{t吨} d日 t吨 + ϕ \sqrt{2 d日 {W公司}_{t吨}} 。

(4.7)

联合分配(4.5)然后修改为

P（P） (θ, λ, ϕ, B类, {x个}_{t吨} ∣ t吨, τ) \propto η^{'} (θ, λ, ϕ) L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） {B类 (t吨)} P（P） ({x个}_{t吨} ∣ λ, ϕ^{2} / λ),

(4.8)

哪里η^′是上的优先分布(θ,λ,ϕ).

为了从联合后验分布中获得蒙特卡罗样本，我们从初始配置开始迭代以下条件采样步骤：

\begin{matrix} θ \sim [θ ∣ λ, ϕ, B类, {x个}_{t吨}, t吨, τ] \propto η^{'} (θ, λ, ϕ) L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}), \\ \begin{array}{l} (λ, ϕ) \sim [λ, ϕ ∣ θ, B类, {x个}_{t吨}, t吨, τ] ⑪ η^{'} (θ, λ, ϕ) P（P） ({x个}_{t吨} ∣ λ, ϕ^{2} / λ), \\ \begin{matrix} B类 \sim [B类 ∣ θ, λ, ϕ, {x个}_{t吨}, t吨, τ] \propto L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） (B类), \\ x个 (t吨) \sim [{x个}_{t吨} ∣ θ, λ, ϕ, B类, t吨, τ] \propto L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） ({x个}_{t吨} ∣ λ, ϕ^{2} / λ) 。 \end{matrix} \end{array} \end{matrix}

The sampling ofθ和(λ,ϕ)是通过MH算法实现的。的条件抽样B类(t吨)和x个_t吨不能以封闭形式实现。我们使用第节中描述的相同的前向算法三以组件方式递归更新它们。

4.3. 缩放变换更新

为了进一步提高计算效率，我们引入了一个更新(λ,ϕ)和x个_t吨同时。这一举措很重要，因为(λ,ϕ)和x个_t吨-从中可以看到方程式（4.7）条件是x个_t吨，后牵引(λ,ϕ)往往较大，而一对较大的(λ,ϕ)反过来，往往会产生较大的x个_t吨因此，如果我们能够采取行动(λ,ϕ)和x个_t吨加在一起，采样器应该能够更快地收敛。

给定的当前配置(λ,ϕ,θ,B类(t吨),x个_t吨)，规模转型提出举措

(ϕ, {x个}_{t吨}) \to (秒 ϕ, 秒 {x个}_{t吨}),

（4.9）

哪里秒是标量。换句话说，该提案试图扩大规模ϕ和控制过程x个_t吨一起向上或向下。保持联合分布(4.8)，标量秒必须根据Liu和Sabatti中的“广义吉布斯抽样”规则从以下条件分布中抽样(2000):

\begin{array}{l} 第页 (秒) \propto 秒^{n个 + 1} P（P） (θ, λ, 秒 ϕ, {B类}_{x个}, {B类}_{年}, {B类}_{z（z）}, 秒 {x个}_{t吨} ∣ t吨, τ) \\ \begin{array}{l} \propto 秒^{n个 + 1} η^{'} (θ, λ, 秒 ϕ) L（左） (t吨, τ ∣ θ, α_{t吨}, 秒 {x个}_{t吨}) P（P） {秒 {x个}_{t吨} ∣ λ, {(秒 ϕ)}^{2} / λ} \\ \propto η^{'} (θ, λ, 秒 ϕ) L（左） (t吨, τ ∣ θ, α_{t吨}, 秒 {x个}_{t吨}), \end{array} \end{array}

(4.10)

其推导被推迟到附录A.由于直接模拟来自配送(4.10)尚不可行，我们应用MH算法时注意到秒=1表示无缩放：首先，建议秒从γ密度Γ(秒;1/c（c）,c（c）)，平均值为1；然后接受提议秒具有类MH概率

第页 = 最小值 {1, \frac{Γ (秒^{- 1}; 1 / c（c）, c（c）) 第页 (秒) 秒}{Γ (秒; 1 / c（c）, c（c）) 第页 (1)}} 。

(4.11)

如果提案被接受，我们将更新(ϕ,x个_t吨)按刻度移动（4.9）。图6比较有刻度移动和无刻度移动的样本的自相关。相对于标准方法的改进是相当大的。从理论角度来看，这种MH移动可以用一种更通用的形式来表述，它可能在其他应用中有潜在的用途。以下结果给出了一般变换群设置下的直接MH规则。证据推迟到附录A。

有无尺度变换更新后验样本的自相关：（a）无尺度更新的λ序列；（b） -无刻度更新的系列；（c） λ-带尺度更新的级数；（d） -带刻度更新的系列

图6

具有和不具有尺度变换更新的后验样本的自相关：（a）λ-无尺度更新的序列；（b）ϕ-无尺度更新的序列；（c）λ-具有比例更新的系列；（d）ϕ-带刻度更新的系列

新标签中打开下载幻灯片

定理1。设∧是状态空间上的局部紧变换群 $X（X）$ ⁠.让μ（d）ν)是∧上的左Haar度量，让π(x个)是感兴趣的目标概率分布 $X（X）$ ⁠假设我们提出以下小组行动：抽签ν从（f）(ν)μ（d）ν)和更新x个^′=ν(x）然后，为了保持π(x个)作为移动后的平稳分布，建议的移动应以概率接受

第页 = 最小值 {1, \frac{（f） (ν^{- 1}) π ({x个}^{'}) | {J型}_{ν} (x个) |}{（f） (ν) π (x个)}},

哪里J型(x个)是转换的雅可比矩阵x个^′=ν(x）。

5.实证研究

5.1. 模拟数据集

为了说明和评估我们的方法的有效性，我们进行了两状态模拟实验。首先，我们生成布朗扩散过程(B类_x个,B类_年,B类_z（z）)，它提供α_t吨根据方程式（3.1）; 然后我们生成两状态马尔可夫过程γ(t吨); 最后是配对 ${({t吨}_{我}, τ_{我})}_{我 = 0}^{n个}$ 根据Cox过程生成（3.2）和指数分布(3.3).使用相同的参数θ，此数据生成过程重复50次，提供50个独立的数据集（如实际实验中的典型情况）。数据集的长度从几百对到几千对不等。为了使模拟数据尽可能真实，模拟中同时包含了背景光子和时间包裹。我们的目标是推断θ从生成的(t吨_我,τ_我).

很容易应用第节中描述的数据增强方法3.2只需将每个实验的似然函数相乘，即可联合分析这50个数据集。对于每个数据集，都会增加一个布朗扩散链。前面有一个平面θ从关节分布中提取5000个后部样本(3.6)。图7显示参数的样本后验分布（垂直条是用于生成数据的真实值）。该算法可以正确识别所有参数。图8成对绘制后验样本。很明显，不需要进一步重新参数化。图9显示了参数的连续蒙特卡罗样本之间的自相关。自相关的快速衰减表明该算法具有快速收敛性。作为进一步的测试，我们从不同的起点启动算法；它们都迅速收敛到后模态区。

布朗扩散（，真值）模拟数据的后验样本直方图：（a）a；（b） b；（c） π1；（d） k；（e） A0型

图7

布朗扩散模拟数据的后验样本直方图(，真值）：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀

新标签中打开下载幻灯片

图8

模拟数据的后验样本配对图

新标签中打开下载幻灯片

模拟数据后验样本的自相关图：（a）a；（b） b；（c） π1；（d） k；（e） A0型

图9

模拟数据后验样本的自相关图：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀

新标签中打开下载幻灯片

为了观察数据增强方法揭示隐藏信息的效果，我们将真实的布朗运动路径与增强路径进行了比较。图10提供了四个典型数据集的详细信息：光曲线（面积）是（多重）增强的布朗因子α(t吨); 真实的α(t吨)如粗曲线所示。数据增强技术似乎很好地恢复了隐藏因素。

模拟数据的增广布朗运动因子（）与实际因子（）的比较：（a）轨迹5；（b）轨迹25；（c）轨迹35；（d）轨迹45

图10

增广布朗运动因子的比较()实际因素()对于模拟数据：（a）轨迹5；（b）轨迹25；（c）轨迹35；（d）轨迹45

新标签中打开下载幻灯片

5.2. 用两态模型分析实际数据集

使用第节中开发的方法三，我们试图分析哈佛大学谢实验室获得的数据集。数据集有1813个观测对。我们获得了5个关键参数的5000个后验样本，其直方图如所示图11在分析中，我们使用信息先验来整合领域（生物化学）知识。例如，从以前的实验中得知，信噪比（即施主光子与背景光子）大约在3/1到8/1之间，这导致我们将宽伽马放在前面A类₀平均值35000，标准偏差25000；拥有k个这远远大于50000实际上不太可能导致指数先验k个平均40000；DNA发夹在打开和关闭状态之间切换的频率大致相等，并且这两个水平一和b条实际上不太可能在之前超过4条线索进入测试版π₁平均值为0.5，标准偏差为0.3，宽伽马先验值为一和b条平均值为2和1，标准偏差为2和0.8。

1813对观测数据集的后验直方图：（a）a；（b） b；（c） π1；（d） k；（e） A0型

图11

1813对观测数据集的后验直方图：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀

新标签中打开下载幻灯片

图12显示了1的后向分布/k个表示开态和闭态之间能量屏障的高度，称为衰减时间常数2.5%和97.5%[39,91]（微秒）为衰减时间常数1提供了95%的对称区间/k个。

图12

衰减时间常数1/k个实验数据集的

新标签中打开下载幻灯片

由于我们的方法使用了似然法，因此它比化学文献中使用的可用的矩型估计方法更有效。在以前的方法中，必须首先将到达时间“装箱”在一起，以消除分子布朗扩散的影响，然后使用装箱到达时间拟合某些矩方程来估计感兴趣的参数（Pfluegl等。,1998; Brown和Silbey，1998). 由于一旦到达时间组合在一起，装箱时间窗口内发生的事情就会丢失，因此装箱方法会严重损失时间分辨率。（从某种意义上说，分箱方法就像是使用某个单位测量距离；如果实际距离比最小单位短，测量将不会给出有意义的结果。）对于我们分析的同一数据集，基于分箱的方法的最大时间分辨率为280μs，并指示衰减时间常数1/k个小于280μs.虽然在质量上与binning方法的结果一致，但我们的分析提供了一个更清晰的推断，即后中位数为1/k个59岁μs和95%对称概率区间[39,91]μ第条。

5.3. 拟合连续扩散模型

由于平稳分布的控制过程(4.2)是N个(0,ξ)，我们注意到两态模型对应于带√的连续扩散模型的退化情况ξ= 0. 这表明，如果后验概率质量√ξ距离0足够远，则数据支持扩散模型。

作为示例，我们应用该算法分析第节中使用的相同模拟数据集5.1，这是根据两国模型生成的。从后向分布中提取了5000个样本(4.8)在1500次迭代的老化周期之后。如所示图13，算法正确识别参数θ和以前一样。此外，√的后验分布ξ表明没有证据表明ξ=0，表示两国模式确实足够。

模拟的两态数据集（，真值）的后验直方图——√ξ的分布高度集中在0附近，强烈表明了两态模型的充分性：（a）a；（b） b；（c） π1；（d） k；（e） A0；（f） √ξ

图13

模拟两状态数据集的后验直方图(，真值）-√的分布ξ高度集中在0左右，强烈表明了两国模式的充分性：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀; （f） √（√）ξ

新标签中打开下载幻灯片

接下来，我们将该方法应用于第节中的实际实验数据集5.2（共有1813个观测对）。比较图14具有图11，我们观察到一,b条和π₁几乎保持不变，而k个和A类₀稍稍移动。转向√的分布ξ，我们注意到模式约为0.25。这种比较似乎提供了一些支持扩散模型的证据，但并不令人信服。因此，我们分析了49个额外的数据集，这些数据集是从对同一DNA发夹的重复单分子实验中获得的，以获得更好的模型评估。

用扩散模型分析实验数据集（1813个数据对）：（a）a；（b） b；（c） π1；（d） k；（e） A0；（f） √ξ

图14

用扩散模型分析实验数据集（1813个数据对）：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀; （f） √（√）ξ

新标签中打开下载幻灯片

图15，显示了八个典型数据集的结果，代表了对单个实验的分析。的结果θ在实验中表现一致，而√ξ对于某些数据集和其他数据集，似乎与0不同。由于所有50个（独立）实验都应由同一组参数控制，我们接下来将联合分析这些数据集（对于每个数据集，布朗扩散和Ornstein–Uhlenbeck控制过程都得到了增强）。图16显示了50个数据集的组合结果。首先，我们注意到，与单个数据集的结果相比，随着更多数据提供更多信息，参数的后验分布变得更窄。第二，更重要的是，图16表明√的后验样本ξ集中在远离0的地方，这强烈表明两国模型不符合数据。

八组实验数据分析结果的方框图：（a）a；（b） b；（c） π1；（d） k；（e） A0；（f） √ξ

图15

八组实验数据分析结果的方框图：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀; （f） √（√）ξ

新标签中打开下载幻灯片

将50个实验数据集（从后验分布中提取5000个样本）一起分析的结果：（a）a；（b） b；（c） π1；（d） k；（e） A0；（f） λ；（g） √ξ

图16

对50个实验数据集进行分析的结果（5000个样本来自后验分布）：（a）一; （b）b条; （c）π₁; （d）k个; （e）A类₀; （f）λ; （g） √（√）ξ

新标签中打开下载幻灯片

虽然每个单独的实验本身只提供了零碎的信息，但将50个实验合并起来表明，两个国家的模型不足以解释实验数据。为了描述DNA发夹分子的详细构象动力学，似乎需要更复杂的模型。

5.4. 贝叶斯因子的蒙特卡罗估计

为了加强图16，我们进一步计算了两国模型和扩散模型之间的贝叶斯因子，这通常用于帮助贝叶斯模型选择（Kass和Raftery，1995). 在我们的例子中，贝叶斯因子可以写成

高炉 = \frac{P（P） (t吨, τ ∣ {M（M）}_{1})}{P（P） (t吨, τ ∣ {M（M）}_{2})} = \frac{\int η (θ) L（左） (t吨, τ ∣ θ, α_{t吨}) P（P） {B类 (t吨)} d日 θ d日 B类}{¦Β η^{'} (θ, λ, ϕ) L（左） (t吨, τ ∣ θ, α_{t吨}, {x个}_{t吨}) P（P） {B类 (t吨)} P（P） ({x个}_{t吨} ∣ λ, ϕ^{2} / λ) d日 θ d日 λ d日 ϕ d日 B类 d日 x个},

哪里M（M）₁表示两国模式M（M）₂表示扩散模型。为了简化表示法，我们将数据表示为年=(t吨,τ)并写入μ=(θ,B类)和ζ=(λ,ϕ,x个_t吨)，使用它可以将上述贝叶斯因子表达式重写为

高炉 = \frac{\int P（P） (年 ∣ {M（M）}_{1}, μ) P（P） (μ ∣ {M（M）}_{1}) d日 μ}{\int P（P） (年 ∣ {M（M）}_{2}, μ, ζ) P（P） (μ, ζ ∣ {M（M）}_{2}) d日 μ d日 ζ},

我们在模型中指出M（M）₂还有一组额外的参数(ζ)与…相比M（M）₁。

由于直接计算这些积分是不可行的，我们寻求通过蒙特卡罗抽样近似贝叶斯因子。我们注意到，如果两个模型的先验分布是一致的，即。

P（P） (μ ∣ {M（M）}_{1}) = \int P（P） (μ, ζ ∣ {M（M）}_{2}) d日 ζ,

这在我们的例子中是正确的，那么贝叶斯因子可以重新表示为似然比的后验均值：

高炉 = \int \frac{P（P） (年 ∣ {M（M）}_{1}, μ)}{P（P） (年 ∣ {M（M）}_{2}, μ, ζ)} P（P） (μ, ζ ∣ 年, {M（M）}_{2}) d日 μ d日 ζ = E类 [\frac{P（P） (年 ∣ {M（M）}_{1}, μ)}{P（P） (年 ∣ {M（M）}_{2}, μ, ζ)} ∣ 年, {M（M）}_{2}] 。

(5.1)

这意味着，如果我们有参数的后验样本(μ^(我),ζ^(我)),我= 1,…,N个，取自较大的模型M（M）₂，我们可以通过以下公式来估计贝叶斯因子

\hat{高炉} = \frac{1}{N个} \sum_{我 = 1}^{N个} \frac{P（P） (年 ∣ μ^{(我)})}{P（P） (年 ∣ μ^{(我)}, ζ^{(我)})} 。

尽管这种方法可能不是最有效的，但它对我们所有的数据都足够有效。

对于第节中分析的特定实验数据集（1813个数据对）5.2，根据5000个后验样本计算得出 $\hat{高炉} = 0.81 \pm 0.10$ ⁠这似乎稍微有利于扩散模型。为了综合可用信息，我们在对50个实验数据集进行组合分析的基础上计算贝叶斯因子，如第节所示5.3，它提供 $\hat{高炉} = (3.43 \pm 0.29) \times 10^{- 9}$ （5000个后验样本的最大似然比为5.79×10⁻⁷). 作为比较，我们还计算了第节中两态模型模拟的50个数据集的贝叶斯因子5.1，我们获得 $\hat{高炉} = 2.51 \pm 0.17$ ⁠这与来自真实实验数据集的贝叶斯因子形成了鲜明对比，这证实了我们对图16。

本小节和上一小节的结果表明，扩散模型更准确地描述了实验数据。从科学的角度来看，这意味着DNA发夹的闭合和开放状态之间的能量屏障比两态模型中描述的简单静态图像具有更复杂的行为。在这种情况下，能量屏障的波动可能是由于DNA分子的构象灵活性，这有待进一步研究。

6.讨论

本文研究了一种单分子实验：DNA发夹分子的荧光寿命涨落。为了有效地分析实验数据，我们讨论了三个重要的统计问题：

（a）
合成和比较用于单分子构象动力学的各种随机模型；
（b）
推导与所涉及的随机模型相关的似然函数；
（c）
解决实验复杂性，尤其是那些涉及未观察到的随机过程的复杂性。

第一个问题是单分子实验本质固有的，与传统的集成实验不同，单分子实验的随机性从根本上需要随机建模。一旦建立了随机模型，有效的推理就依赖于似然公式；我们演示了如何使用离散化和矩阵技术来获得闭合形式的可能性。分析真实实验数据的一个重要方面是处理实验并发症。一些复杂情况（例如时间包装和延迟时间的负读数）相对容易处理，对可能性的修改也很简单。由于无法获得可计算的分析表达式，其他复杂情况，如与未观察到的随机过程相关的复杂情况，带来了更大的挑战。事实证明，借助马尔可夫链蒙特卡罗方法的数据增强技术是一种非常强大的工具。它们不仅概念简单，而且提供了一种可行的方法来规避分析的棘手问题。

我们对贝叶斯因子的蒙特卡罗估计依赖于简单直观的恒等式(5.1)似然比的后验均值（在较大模型下）等于两个嵌套模型的贝叶斯因子。这个恒等式似乎是众所周知的（我们的两位同事，孟晓丽和唐纳德·鲁宾，以前曾私下使用过这种恒等式的变体），并且与用于推导桥式抽样的公式密切相关（孟和王，1996)以及比率重要性抽样方法（Chen和Shao，1997; 陈等.,2000). 然而，我们既没有找到明确说明这个结果的引用，也没有看到它在蒙特卡罗计算中的大量使用。在较小模型为固定分布（无未知参数）的特殊情况下，我们的方法等效于牛顿和拉弗瑞的“调和平均”方法(1994). 尽管我们的方法在特定情况下可能不如其他一些更精细的方法有效（例如，陈中描述的技术等。(2000)孟和王(1996))，它对本文所报道的研究有效——在所有情况下，似然比的方差都得到了很好的控制。相比之下，Chib的(1995)该方法在这里不起作用，因为它需要估计某些点的后验密度值，这对于我们的模型是不可行的。可逆跳跃法（Green，1995)也不太可能有帮助，因为我们知道没有有效的跳跃规则来添加或删除隐藏的Ornstein–Uhlenbeck过程（即在两个状态和连续扩散模型之间跳跃）。

尽管本文所述的马尔可夫链蒙特卡罗方法在过去十年中已在统计界得到广泛接受（刘，2001)，它们在其他科学学科中的使用相对较少。在过去，许多物理科学研究人员认为没有必要使用精细而有效的统计推断方法，因为他们在集合上的数据通常非常庞大特别的诸如矩匹配之类的方法将足以为他们提供所需的信息。现代技术的进步所促成的单分子实验，以及许多大规模基因组学实验，似乎大大改变了环境。如本文所示，与现有的矩匹配和装箱方法相比，贝叶斯分析对与DNA发夹动力学相关的参数提供了更清晰的估计。数据增强方法还使我们能够研究不同模型与实验数据的拟合优度，这对于现有的基于矩的方法来说尤其困难，甚至是不可能的。即使在二十年前，单分子实验和数据增强（使用马尔可夫链蒙特卡罗）方法都是不可想象的。因此，我们希望本文能进一步激发人们对将现代统计方法应用于有趣和重要的科学问题的兴趣。

本文中使用的实验数据集可以从作者的网页上下载网址：http://www.fas.harvard.edu/~skou（斯库）和网址：http://www.fas.Harvard.Edu/~俊柳。

寇、谢和刘的论文讨论

阿兰·霍克斯(威尔士大学，斯旺西)

我欢迎随机模型论文在学会的一次会议上出现，虽然重点可能更多地放在问题的推断方面，但这种出现可能太罕见了。对于没有扩散的简单模型，我倾向于使用好的老式似然法。然而，一旦扩散进入模型，我就能很好地理解使用马尔可夫链蒙特卡罗方法的贝叶斯方法的吸引力：如果你要用随机生成的值来增加数据，那么你也可以生成参数值。我对马尔可夫链蒙特卡罗抽样的实用性不够精通，无法对其细节作出认真的评论，尽管它似乎取得了令人印象深刻的效果。从现在开始，我将仅限于讨论随机建模。

我想附和作者关于单分子观测在理解化学过程的结构和功能方面的优势的介绍性发言。我与伦敦大学学院的David Colquhoun和其他人一起参与了单离子通道行为的建模大约25年。我想简单地提请您注意。霍克斯作了初步介绍(2003).

离子通道是嵌入细胞膜中的大型蛋白质分子，控制带电离子在膜上的运动。这些带电离子的流动构成了电流的流动。神经系统中的所有电活动似乎都受到离子通道的调节，因此在许多不同的活动中发挥作用，包括思维过程、神经信号的传递及其转化为肌肉收缩和控制胰岛素的释放。了解他们的行为有助于我们了解正常生理以及药物和毒素的作用。因此，这是朝着开发各种医疗条件的治疗方法迈出的重要一步。

值得注意的是，自从尼赫和萨克曼获得诺贝尔奖以来(1976)，可以观察到流过单个通道的电流。在许多情况下，除了一点噪声外，当通道关闭时，电流为零，或者当通道打开时，电流是恒定值。提取一系列交替开闭周期的理想信号通常并不太困难，如图17实际上，我们很幸运能够比荧光实验中更直接地观察到这个分子。

图17

从噪声信号中提取的理想开-闭序列

新标签中打开下载幻灯片

与本文一样，基本模型是一个有限状态马尔可夫过程，它描述了分子可以采用的各种物理状态。每个状态可分为打开或关闭。图18显示了甘氨酸受体的最新可能模型（Burzomato等。,2004)有七个关闭状态和三个打开状态（图18). 从任何这样的模型中，我们都可以用（隐藏的）半马尔可夫过程来描述开-闭序列。我们观察到的是一个聚集的马尔可夫过程，因为我们只能观察它处于哪一组状态（打开或关闭），而不能观察单个状态。实际上，有一些关于区间太短而无法解决的复杂问题：半马尔可夫公式可以考虑到这些问题，但我们不需要在这里详细介绍。由此，我们可以推导出这样一个序列的似然性，并且它已经被证明在获得参数的良好最大似然估计和区分潜在的潜在模型方面是非常有效的。

图18

甘氨酸受体的可能模型：最下面一行的三个状态是开放状态；其他七个州关闭

新标签中打开下载幻灯片

贝叶斯方法也被用于直接从原始信号加噪声测量中得出参数估计值，以及从噪声中提取信号（Fredkin和Rice，1992). 即使在多层次观测和高噪声水平下，这些方法也能很好地工作（参见中的示例图19).

图19

贝叶斯恢复：（a）四级模拟信号；（b）信号加噪声；（c）恢复信号（四个非常短的事件丢失）

新标签中打开下载幻灯片

在转移注意力之后，让我回到本文中的模型。

我对基本的随机模型感到困惑，毫无疑问，这是因为我对物理的掌握不够。

如果你看一下可能性的推导方程式（2.3）你看，光子观测的过程取决于潜在的伽马过程γ(t吨)它本身并不以任何方式依赖于脉冲序列。然而，我们被告知是脉冲刺激光子发射：令人困惑！

然后，延迟时间似乎在概率上独立于光子到达过程，并且彼此独立，并且在数学上取决于到达时伽马过程的值，t吨_我仅限。乍一看，这些假设似乎不太可能。我想答案是，如果脉冲速率远快于γ(t吨)变化和光子到达的速率。

我本希望能对这些问题进行更多的解释。

有几个小问题。

在第474页的备注1中，我们被告知推导可能性的离散化方法，如附录A与从一组包含无穷小生成器的微分方程中导出可能性相比，更容易推广到更复杂的潜在马尔可夫模型和布朗扩散的引入。我不同意。在我看来，微分法更简单、更优雅。此外，使用这种方法，更复杂的模型实际上并不比简单模型更难求解。下面是沿着这些线的推导。

稍微概括一下，让我们X（X）(t吨)表示分子当时的状态t吨并假设此时的泊松强度为α(t吨)λ_我如果X（X）(t吨) =我。

让H（H）=诊断(λ₁,λ₂，…，λ_k个)然后让矩阵G公司(t吨)具有元素

克_{我 j个} (t吨) = P（P） {Y（Y） (t吨) - Y（Y） (0) = 0 \cap X（X） (t吨) = j个 ∣ X（X） (0) = 我} 。

那么很容易看出

\begin{array}{l} G公司 (t吨 + δ t吨) = G公司 (t吨) {我 - α (0) H（H） δ t吨 + o个 (δ t吨)} {我 + 问 δ t吨 + o个 (δ t吨)} \\ = G公司 (t吨) [我 + {问 - α (0) H（H）} δ t吨 + o个 (δ t吨)] \end{array}

以便

\frac{d日 G公司}{d日 t吨} = G公司 {问 - α (0) H（H）},

这很容易解决

G公司 (t吨) = 经验 [{问 - α (0) H（H）} t吨]

根据需要。

我还有一个小小的诡辩。在节中4.1我们被告知，在二乘二模型中A类₁和A类₂实验上无法区分这对B类₁和B类₂从而形成所谓的聚集马尔可夫过程。只有当常量一₁=一₂和b条₁=b条₂在以下定义中方程式（2.4）。

我很高兴提议投票感谢作者们提交的激励性论文。

Mark钢材(考文垂华威大学)

值得祝贺的是，作者们发表了一篇非常有趣且富有启发性的论文，该论文建立在有关单分子实验的生物化学新兴文献的基础上。恐怕我在生物信息学方面没有真正的专业知识，因此，我将不得不将此讨论局限于贝叶斯统计中更为标准的问题。

具有潜在过程的贝叶斯推理

作者提出的模型依赖于未观察到的随机过程，这无疑会使计算分析复杂化。关于贝叶斯推断与潜在随机过程的文献越来越多（例如在金融领域），我认为，将这里使用的方法与文献中的一些备选方案进行比较是很有意思的。

本文的第一个问题是作者使用了扩散链的组件式更新。在第节中，提到了“阻塞”扩散的采样器3.2以及在后面的章节中回到这些分块移动的承诺，但这可能是指团队移动，而不是“真正的”分块移动，例如刘等。(1994)或者Shephard和Pitt(1997). 在这种情况下，更常见的阻塞类型是否有帮助？

第二个可能更关键的问题是如何处理过度调节（即潜在变量和参数之间的高度相关性）。这种情况发生在x个(t吨)在连续扩散模型中引入的过程。提出的解决方案是使用一个尺度变换，该变换联合更新过程及其扩散参数。这种规模转换可能存在的一个问题是，整个过程是一次性缩放的，如果这种可能性对过程的某些部分特别有用，那么这种移动可能很难接受。在这里，它似乎奏效了，可能是刘和萨巴蒂的“广义吉布斯”步骤的结果(2000)，这就发挥了可能性。这种缩放有点像罗伯茨和斯特拉默的方法(2001)世卫组织对部分观测扩散模型中的“缺失数据”进行了重新参数化（尽管他们使用了稍微不同的转换来实现观测点处样本路径的连续性）。这可以解释为所谓的“非中心参数化”，如Papaspiliopoulos所述等。(2003). 基本上，非定中心解决了增强数据包含的参数信息比观测数据多得多的情况。巴恩多夫·尼尔森和谢泼德(2001)介绍了一类随机波动率模型，其中（潜在）波动率被建模为Ornstein–Uhlenbeck过程，由无高斯分量的正Lévy过程驱动。在这些模型的背景下，还提出了非中心化来解决Roberts中的过调节问题等。(2004). 对于相同类型的模型，Griffin和Steel(2003)提出了一种根据提出的参数值对波动过程进行稀释的采样器。格里芬钢铁公司(2003)实现了所谓的“依赖细化”，即在过程中添加或删除的跳跃往往相对较小。这是通过使用Lévy跳跃过程的Ferguson–Class表示实现的，它允许我们在提出过程更改时关注小跳跃。因此，新流程将与前一流程“接近”，同时也与参数的新值兼容。事实上，罗伯茨等。(2004)提到这种方法也可以解释为非中心参数化。尽管上述论文中使用的模型有些不同，但我想知道在本论文中是否可以使用类似的想法。

以前的问题

作者在第节中提到1“我们对涉及的各种参数有相当详细的了解（事先）”。然而，本文很少讨论先验启发或一般先验问题。在模拟数据示例中，始终使用平面先验。作者没有明确“平坦”是什么意思，但如果它意味着统一（对数）在ℜ上定义的参数₊这使得先验不恰当，并引发了人们熟悉的后验存在的幽灵。后验与所选的前验存在吗？此外，为什么使用的特定参数化的一致性（而不是其他参数化）是反映缺乏先验信息的合理选择（这大概是作者的意图）？是否有一个不变性论证可以激发这种选择？此外，通过例如事先敏感性分析，了解事先假设的重要性也很有用。虽然典型大小的数据集可能会在某些方向上淹没先验数据，但我不认为所有先验维度都必然如此。最好知道在引出实际问题的先验时，应该将精力集中在哪里，而关于如何进行此类引出的一些想法对于这些方法的应用用户来说可能非常有价值。

贝叶斯因子的估计

中的身份方程式（5.1）用于估计贝叶斯因子的方法很有意思，正如作者所说，在文献中没有太多使用。我想知道原因是否与“调和平均值”估计值相同，如果更简单的模型没有未知参数，则构成特殊情况。牛顿和拉弗瑞(1994)评论说，调和平均值通常不满足高斯中心极限定理，因此可能非常不稳定。我们预计，对于似然比来说，这样的问题会稍微不那么严重，但我想知道由此产生的估计的稳定性是否仍然是一个问题，特别是如果更大模型中的额外参数涉及潜在过程。

最后，我很高兴再次投票表示感谢。

感谢票以鼓掌方式通过。

克里斯·夏洛克(兰卡斯特大学)

我们评论了作者在他们的模拟研究中使用平坦先验信息，并提供了吉布斯采样器（参见费恩黑德和夏洛克(2004))作为第节中描述的大都会黑斯廷斯计划的替代方案2。

平坦先验

我们观察到，在某些参数的整个范围内，可能性的下限严格大于0。例如

P（P） (数据) ⩾ P（P） (数据 \cap 分子总是闭合的)

其中右侧独立于参数b条。

因此，（在标准或对数刻度上）b条会导致后路不合适。

吉布斯采样器

在当前背景下隐（马尔可夫）链的状态DNA分子是闭合的还是开放的。

吉布斯采样器有三个阶段。

第1步：在上一次迭代结束时给定参数值，从其精确条件分布中采样光子发射时间的隐马尔可夫过程状态。

第2步：对于排放之间的每个时间间隔，给定链在端点的状态，从这些链的精确条件分布中对链进行采样。

步骤3：给出完整潜在马尔可夫链的这个实例，从给定数据的精确参数分布中采样。

步骤1使用向前-向后算法，如果我们使用共轭先验，步骤3很简单。

第二步是基于费恩黑德和梅利格科茨杜的想法(2004)并涉及到创建一个支配泊松过程，其速率是状态变化过程和光子发射过程强度之和的最大值（链的所有状态）。这些主要事件对应（通过变薄）可能的状态变化或光子发射。在两个光子发射之间的每个间隔内，我们模拟

（a）
根据其真实条件分布得出的主导事件数，
（b）
每个主导事件的时间，在区间内独立且一致，以及
（c）
每个主要事件的状态变化（如有）；这里我们再次使用了前向-后向算法。

我们重复了作者在章节中的模拟2（带有模糊但适当的先兆）。图20显示了两个参数的结果。参数k个₂₁在五个参数中表现出最差的混合。模拟所花费的中央处理器单元时间不到带有块更新的加性随机遍历Metropolis方案的两倍。

（a），（b）迹图，（c），（d）密度估计（，真实参数值）和（e），（f）（a），（c），（e）λ1（：=A0/a）和（b），（d），（f）k21的自相关函数

图20

（a），（b）迹线图，（c），（d）密度估算(，真参数值）和（e），（f）（a），（c），（e）的自相关函数λ₁(: =A类₀/一)和（b）、（d）、（f）k个₂₁

新标签中打开下载幻灯片

奥米洛斯·帕帕斯皮利奥普洛斯(兰卡斯特大学)

值得祝贺的是，作者发表了这篇有趣的论文。我想对计算方法作一些评论，并通过给出本文中缺少的一些关键参考文献，提请作者注意隐随机过程推理方面的重要工作。

第节中描述的模型2是一个马尔可夫调制泊松过程，在拟合这些模型方面已经做了很多工作；参见Asmussen示例(2000)以及其中的参考。很有可能方程式（2.3）是标准的，可以派生（以及方程式（3.11）)使用比中采用的方法简单得多的（矩阵分析）参数附录A。

如果w个_xy公司=w个_z（z）在里面方程式（3.1）（或者如果我们愿意假设这个过程存在于两个维度中）我们可以避免增加三个布朗运动，而直接增加α(t吨)该过程是已知跃迁密度的贝塞尔过程的变换。这可能导致更高效的马尔可夫链蒙特卡罗采样，因为有大量关于α，而这三个过程是部分不可识别的（所有距离原点相同的路径都具有相同的可能性）。作者是否因路径空间中的多种模式而遇到混合问题？第节中建议的隐藏时间序列的逐时间点更新3.2虽然它会导致快速更新，但众所周知，它会导致非常缓慢的混合马尔可夫链算法。皮特和谢泼德(1999)基于罗伯茨和萨胡的研究结果(1997)，给出了一些解析收敛结果；另请参阅提花机中的讨论等。(1994)随机波动率模型的经验结果，例如Elerian等。(2001)关于更新扩散路径的有效阻塞方案的建议。在如此高的维度上，很难用图形来诊断混合效果是否良好，但我仍然不太相信作者的观点图10，因为在许多情况下，“真实”路径位于模拟路径覆盖的范围之外。

罗伯茨和斯特雷默(2001)表明吉布斯采样器更新了扩散数据X（X）_t吨=b条(X（X）_t吨)d日t吨+σ(X（X）_t吨,ϕ)d日W公司_t吨以及参数θ收敛性差；事实上，如果所有路径都增加了，那么采样器是可以减少的。他们建议进行非中心化的重新参数化

\begin{matrix} {X（X）}_{t吨} \to {\tilde{X（X）}}_{t吨} = η ({X（X）}_{t吨}; ϕ), \\ η (x个; ϕ) = \int_{0}^{x个} \frac{1}{σ (u个; ϕ)} d日 u个; \end{matrix}

另见罗伯茨等。(2004)，帕帕斯皮利奥普洛斯等。(2003)和帕帕斯皮利奥普洛斯(2003). 重新参数化和边缘增强方法之间有着密切的联系，如第节中使用的方法4.3在这种情况下，前一种方法更加灵活，因为它可以很容易地处理任意扩散函数和离散观测的扩散，但了解这两种方法是否可以组合是很有趣的。

马丁·雅各布森和迈克尔·瑟伦森(哥本哈根大学)

我们祝贺作者发表了有趣且富有启发性的论文，并希望指出在更一般的情况下（没有任何额外的复杂性）对其似然函数进行了严格推导。

假设在时间间隔[0，T型]. 数据是成对的事件时间和标记(t吨₁,τ₁),…,(t吨_n个,τ_n个) (0<t吨₁<…<t吨_n个<T型). 计数过程的强度为λ_t吨=λ(M（M）_t吨−)，其中M（M）是状态空间为{1，…，的连续时间马尔可夫链，…，d日}，强度矩阵问和平稳分布π₁,…,π_d日标记在某个时间的概率密度t吨以标记计数过程的过去和马尔可夫链为条件M（M）是（f）(τ|t吨,M（M）_t吨−). 马尔可夫链上有标记计数过程的似然函数M（M）是

{L（左）}_{c（c）} = {\prod_{我 = 1}^{n个} λ_{{t吨}_{我}} （f） (τ_{我} ∣ {t吨}_{我}, {M（M）}_{{t吨}_{我} -})} 经验 (- \int_{0}^{T型} λ_{秒} d日 秒)

（关于适当的主导措施）；参见Bremaud示例(1981). 迭代定义

{L（左）}_{我} (米) = {E类}_{M（M）} {{Z轴}_{我} {L（左）}_{我 + 1} ({M（M）}_{{t吨}_{我}}) ∣ {M（M）}_{{t吨}_{我 - 1}} = 米}, 我 = 1, \dots, n个 + 1,

哪里E类_M（M）表示对以下分布的期望M（M）,L（左）_n个+2(米)=1和

{Z轴}_{我} = {q个}_{我} ({M（M）}_{{t吨}_{我}}) 经验 {- \int_{{t吨}_{我 - 1}}^{{t吨}_{我}} λ ({M（M）}_{秒}) d日 秒},

具有t吨₀= 0,t吨_n个+1=T型,q个_n个+1(米)=1和q个_我(米) =λ(米)（f）(τ_我|t吨_我,米). 然后，通过迭代条件期望并使用马尔可夫特性，我们得到了无条件似然

L（左） = {E类}_{M（M）} ({L（左）}_{c（c）}) = {E类}_{M（M）} ({L（左）}_{1}) = (π_{1}, \dots, π_{d日}) {L（左）}_{1}

哪里L（左）_我=(L（左）_我（1），…，L（左）_我(d日))^T型（T表示换位）。使用第页_我(米) =q个_我(米)L（左）_我+1(米),

\begin{array}{l} {L（左）}_{我} (米) = {E类}_{M（M）} (\underset{K（K） \to \infty}{林} [{第页}_{我} ({M（M）}_{{t吨}_{我}}) 经验 {- \sum_{j个 = 1}^{K（K）} λ ({M（M）}_{秒_{j个}}) 小时} ∣ {M（M）}_{{t吨}_{我 - 1}} = 米]) \\ = \underset{K（K） \to \infty}{林} ({E类}_{M（M）} [{第页}_{我} ({M（M）}_{{t吨}_{我}}) \prod_{j个 = 1}^{K（K）} 经验 {- λ ({M（M）}_{秒_{j个}}) 小时} ∣ {M（M）}_{{t吨}_{我 - 1}} = 米]) \end{array}

哪里秒_j个=t吨_我−1+j小时,小时= (t吨_我−t吨_我−1)/K（K），以及我们在哪里使用了支配收敛定理。使用与作者相同的矩阵计算，我们得到

\begin{matrix} {L（左）}_{我} (米) = \underset{K（K） \to \infty}{林} [\sum_{米_{1}, \dots, 米_{K（K）}} {第页}_{我} (米_{K（K）}) \prod_{j个 = 1}^{K（K）} 经验 {- λ (米_{j个}) 小时} 第页 (米_{j个} ∣ 米_{j个 - 1})] \\ \begin{array}{l} = \underset{K（K） \to \infty}{林} [\sum_{米_{1}, \dots, 米_{K（K）}} {第页}_{我} (米_{K（K）}) \prod_{j个 = 1}^{K（K）} {经验 (问 小时) 经验 (- H（H） 小时)} 米_{j个 - 1} 米_{j个}] \\ = {[经验 {(问 - H（H）) ({t吨}_{我} - {t吨}_{我 - 1})} H（H） {D类}_{我} {L（左）}_{我 + 1}]}_{米}, \end{array} \end{matrix}

对于我≼n个，其中第页(米_j个|米_j个−1)是马尔可夫链的转移概率，米₀=米,H（H）=诊断{λ(1),…,λ(d日)}和D类_我=诊断{（f）(τ_我|t吨_我,1),…,（f）(τ_我|t吨_我,d日)}. 总之，

{L（左）}_{我} = 经验 {(问 - H（H）) ({t吨}_{我} - {t吨}_{我 - 1})} H（H） {D类}_{我} {L（左）}_{我 + 1}

对于我≼n个、和L（左）_n个+1=经验{(问−H（H）)(t吨_我−t吨_我−1)}e（电子），其中e（电子）=(1,…1)^T型，所以通过迭代

L（左） = (π_{1}, \dots, π_{d日}) [\prod_{我 = 1}^{n个} 经验 {(问 - H（H）) ({t吨}_{我} - {t吨}_{我 - 1})} H（H） {D类}_{我}] 经验 {(问 - H（H）) (T型 - {t吨}_{n个})} e（电子） 。

假设函数λ(米)是确定性的，但对证据的检查表明，当允许强度取决于前一事件之前标记计数过程的行为时，它也适用，而在事件之间，它还取决于M（M）同上。强度也可能取决于其他随机性来源，因此本文中的其他模型也包括在内，但重要的是，事件之间对此类额外随机性的依赖保持不变。时间间隔内的观察[t吨₁,t吨_n个]被服用覆盖t吨₁=0和T型=t吨_n个。

理查德·J·男孩(纽卡斯尔大学)

值得祝贺的是，作者发表了一篇令人印象深刻的论文。生物化学实验的最新进展现在提供了对单个细胞内生物机制的洞察。与微阵列和聚合酶链反应等集合平均实验相比，这些技术可能会让我们更全面地了解细胞生物学过程。在纽卡斯尔，达伦·威尔金森（Darren Wilkinson）和我正在开发马尔可夫链蒙特卡罗（Markov chain Monte Carlo）方法，通过使用关于单个细胞内化学物种分子数量的部分观测时间进程数据推断生化网络中的随机动力学常数。这个工作（男孩们等。,2004; 威尔金森等。,2004)还使用了数据增强，在我们的例子中，对由于部分时间过程信息而产生的额外不确定性进行平均。

斯蒂尔教授对模型敏感性的评论，包括对先验信息的使用，很有意思。评估先验信息的影响可能非常困难，特别是当在具有潜在结构的模型中使用扩散先验时，这些潜在结构通常只能通过数据微弱地识别。在更简单的场景中，可以通过观察以下边际可能性来评估对先前规范的敏感性

π (年 ∣ μ) \propto \hat{π} (μ ∣ 年) / π (μ),

哪里 $\hat{π}$ 是根据马尔可夫链蒙特卡罗输出计算的边际后验密度的经验估计。事实上，该分析使用的是扩散先验，模型似乎不存在可识别性问题，因此边缘后验几乎是标准化的边缘似然。谈到模型敏感性，评估用于Ornstein–Uhlenbeck过程的离散化效果的一种方法是，通过插入时间点来考察推论对具有更精细离散化的模型的敏感性t吨_我,1,…,t吨_我,克然后使用数据增强技术整合出附加值x个(t吨_我,j个).

作者评论说，他们之前没有发现贝叶斯因子计算的参考(5.1)这是基于假设似然比的后验平均值H（H）₀以下为：ζ= 0. 艾特金等。(2005)，我们研究如何使用这种后验分布来评估这些假设。不只是使用平均值，而是使用完全分布来衡量与假设相反的证据强度。很有意思的是，在这个分析中，似然比的后验分布是否特别偏斜，以及是否有一些可变性和分位数的度量。

这个作者随后以书面形式回复如下。

我们感谢讨论者对我们的论文发表了深思熟虑的建设性意见。我们将首先讨论最常见的问题，然后对特定讨论者提出的问题做出（必要的）简短回答。

贝叶斯推理、先验规范和模型敏感性

对我们来说，随机模型的贝叶斯推断的一个主要优点是，所有量，无论是观测数据、未观测到的潜在过程还是感兴趣的参数，都在一致的概率框架下得到了相同的基础。因此，贝叶斯方法在概念上直接用于推断具有潜在过程的随机模型，并已应用于金融文献中研究部分观测到的扩散过程、随机波动率模型和期限结构模型，正如Steel和Papaspiliopoulos所指出的那样。尽管所有这些模型都属于马尔可夫或隐马尔可夫模型的一般范畴，但我们的模型类似于更经典的隐马尔可夫模型，而不是部分观察到的扩散过程。

贝叶斯分析的一个组成部分是了解先前规范的影响。在分析大多数自然科学实验时，先验分布应尽可能分散，这样科学家的数据才能得出主要的推断结论。一个“好的”实验通常是指一个提供信息的实验，很容易压倒先前的分布。在我们的分析中，我们对所有参数使用“平坦”先验，其中平坦性仅意味着在有限但相对较大的区域上均匀分布，这由领域科学知识决定。虽然我们的统一先验不是变换不变的，但我们的灵敏度分析表明，先验的影响确实非常小。用于第节中的模拟研究5.1，我们生成了50个数据集，每个数据集包含数百到数千个数据对。对于广泛的先验分布，后验分布基本上没有差异。

Sherlock提出的关于我们的两态模型的似然函数奇异性的优秀观点与典型高斯混合模型中的相似，并且可以通过使用适当的先验在理论和实践上避免。

男孩们提出了一种有趣的方法来测试模型敏感性，特别是对于扩散的Ornstein–Uhlenbeck模型。我们只是在论文中证明了简单的两态模型不足以描述DNA发夹动力学。Boys建议的模型敏感性测试可能有助于进一步测试扩散模型的有效性。

可能性计算

我们感谢雅各布森和瑟伦森对可能性进行了严格推导。对我们来说，离散化和矩阵方法非常直观，因为它提供了一种将处理离散时间过程的技术扩展到连续时间过程的自然方法。离散化方法的另一个好处是，它可以为分析和计算问题提供答案或见解：当似然函数的分析形式存在时，可以通过让离散化长度收缩到零来获得它（如论文所示），而如果没有分析解，该方法对应于计量经济学文献中使用的数值“欧拉方法”，是最近用于推断部分观测扩散过程的数据增强方法的基础（Elerian等。,2001).

马尔可夫链蒙特卡罗问题

有一些关于使用“块”移动来增强潜在过程的建议。我们同意这在一般状态空间模型中非常有用，特别是当我们可以很容易地从整个变量块的联合条件分布中取样时（Carter和Kohn，1996). 当这种条件块分布无法轻松处理时，Shephard和Pitt(1997)提出了一些可以提高计算效率的部分更新方法，如添加移动。但根据我们的经验，与真正的块移动相比，这些部分更新的效果要小得多。在刘和萨巴蒂(2000)，我们为部分观测到的扩散过程引入了一种多级部分块更新方法，它进一步改进了Shephard和Pitt的方法，但改进仍然不够大，不值得为我们当前的问题进行额外的编程工作。

然而，我们同意帕帕斯皮利奥普洛斯的观点，即我们的三维布朗路径是不可识别的——在布朗路径围绕z（z）-轴。因此，我们的单站点更新采样器无法充分探索整个路径空间。我们的方法似乎仍然有效的原因是，可能性仅取决于α(t吨)取样器在角度空间中是否混合良好并不重要。在不变性论证的基础上，我们还可以引入一个新的步骤来加快收敛速度：B类^′(t吨) =对○B类(t吨)，其中对是一种类似旋转的变换α(t吨)保持不变。这样的动作应该容易被接受，并且应该比阻挡动作效果好得多。

我们也同意Steel和Papaspiliopoulos的观点，即重新参数化可以显著改善收敛性。我们觉得组移动可能比大多数重新参数化更直观、更容易使用，因为

（a）
它可以实现几乎所有的重新参数化效果，而无需实际改变原始分布和
（b）
我们可以在一次马尔可夫链蒙特卡罗迭代中轻松地进行不同类型的群移动，但同时进行几种类型的重新参数化要复杂得多。

然而，一些重新参数化，例如Roberts和Stramer(2001)不能简单地表述为群体运动，这表明可能存在一个比广义吉布斯结构更有用、更通用的数学框架（Liu和Sabatti，2000)包含非组转换。

实验背景

我们感谢霍克斯指出当前工作与建模和推断离子通道行为的文献之间的联系，并就光子到达时间和延迟时间提出澄清问题。在我们的实验中，由于激光脉冲速率比γ(t吨)正如霍克斯正确观察到的那样，光子到达率本质上是由分子的潜在状态决定的。激光功率影响常数A类₀光子延迟时间测量分子从被激发时发射光子所需的时间长度；这个时间长度基本上与激光脉冲速率无关，并且由分子的状态决定。在实验中，将分子置于激光焦体积中。焦体积是椭圆形的，在x个-和年-轴，但在z（z）-轴，导致w个_xy公司=310 nm和w个_z（z）=1760 nm英寸方程（3.1）因此，一个贝塞尔过程不足以增加α(t吨)正如帕帕斯皮利奥普洛斯所建议的那样。

最后，我们感谢所有讨论者的深刻贡献。我们希望读者会喜欢阅读这些评论，并思考我们提出的各种科学、统计和计算问题。

致谢

作者感谢Haw Yang提供的实验数据和图1和三以及Long Cai和Xiao-Li Meng进行了有益的讨论。作者感谢裁判，他们的意见和建议帮助我们澄清了许多问题，改进了整体表现。这项工作得到了两项国家科学基金会拨款（DMS-0204674和DBI-0138028）、哈佛大学克拉克库克基金和能源部科学办公室基础能源科学办公室的拨款（DOE-FG02-00ER15072）的部分支持。

工具书类

阿格蒙

,

N。

和

霍普菲尔德

,

J·J。

(

1983

)

垂直于反应坐标的有界扩散化学反应的瞬态动力学：缓慢构象变化的分子内过程

。

化学杂志。物理学。

,

78

,

6947

–

6959

。

谷歌学者

交叉参考

书目数据库

安萨里

,

答：。

,

库兹涅佐夫

,

S.V.公司。

和

沈

,

年。

(

2001

)

折叠漏斗的构型扩散描述了DNA发夹的动力学

。

程序。国家。阿卡德。科学。美国

,

98

,

7771

–

7776

。

谷歌学者

交叉参考

书目数据库

安萨里

,

答：。

,

沈

,

年。

和

库兹涅佐夫

,

S.V.公司。

(

2002

)

错误折叠的环降低了DNA发夹形成的有效率

。

物理学。修订稿。

,

88

,

069801

。

发动机罩

,

G.公司。

,

克里切夫斯基

,

O。

和

利布沙贝

,

答：。

(

1998

)

DNA发夹状突起构象波动动力学

。

程序。国家。阿卡德。科学。美国

,

95

,

8602

–

8606

。

谷歌学者

交叉参考

书目数据库

棕色

,

F.L.公司。

和

西尔比

,

R·J。

(

1998

)

低温玻璃中两能级系统耦合对单分子线型影响的研究

。

化学杂志。物理学。

,

108

,

7434

–

7450

。

谷歌学者

交叉参考

书目数据库

曹

,

J。

(

2000

)

单分子动力学的事件平均测量

。

化学。物理学。莱特。

,

327

,

38

–

44

。

谷歌学者

交叉参考

书目数据库

陈

,

男-女。

和

邵

,

问-米。

(

1997

)

估计归一化常数比值的蒙特卡罗方法

。

安。统计师。

,

25

,

1563

–

1594

。

谷歌学者

OpenURL占位符文本

书目数据库

陈

,

男-女。

,

邵

,

问题-问题。

和

易卜拉欣

,

J·G·。

(

2000

)

贝叶斯计算中的蒙特卡罗方法

。

纽约

以下为：

施普林格

。

芯片

,

美国。

(

1995

)

吉布斯输出的边际似然

。

《美国统计杂志》。助理。

,

90

,

1313

–

1321

。

谷歌学者

交叉参考

书目数据库

考克斯

,

D.R.公司。

(

1955

)

包含补充变量的非马尔可夫随机过程分析

。

程序。外倾角。菲尔·索克。

,

51

,

433

–

441

。

谷歌学者

交叉参考

书目数据库

埃格林

,

C、。

,

炸薯条

,

J。

,

品牌

,

L。

,

Günther公司

,

R。

和

塞德尔

,

C、。

(

1998

)

用选择性荧光光谱法监测单个分子的构象动力学

。

程序。国家。阿卡德。科学。美国

,

95

,

1556

–

1561

。

谷歌学者

交叉参考

书目数据库

埃利奥特

,

R。

,

阿贡

,

L。

和

摩尔

,

J。

(

1997

)

隐马尔可夫模型：估计与控制

。

纽约

以下为：

施普林格

。

弗雷德金

,

D。

和

大米

,

J。

(

1986

)

关于聚集马尔可夫过程

。

J.应用。普罗巴伯。

,

23

,

208

–

214

。

谷歌学者

交叉参考

书目数据库

弗罗利奇·阿蒙德

,

美国。

,

大风

,

英国。

和

Osheroff公司

,

N。

(

1994

)

拓扑异构酶II对DNA发夹的定点切割。DNA二级结构作为酶识别/裂解的决定因素

。

生物学杂志。化学。

,

269

,

7719

–

7725

。

绿色

,

P.J.公司。

(

1995

)

可逆跳跃马尔可夫链蒙特卡罗计算与贝叶斯模型确定

。

生物特征

,

82

,

711

–

732

。

谷歌学者

交叉参考

书目数据库

格伦威尔

,

J。

,

玻璃

,

J。

,

拉科斯特

,

T。

,

德尼

,

答：。

,

尚拉

,

D。

和

舒尔茨

,

第页。

(

2001

)

利用单对荧光能量转移监测DNA发夹的构象波动

。

美国化学杂志。Soc公司。

,

123

,

4295

–

4303

。

贾

,

年。

,

锡特尼克

,

答：。

,

锂

,

L。

,

弗拉迪米洛夫

,

美国。

,

库珀曼

,

B。

和

霍克斯特拉瑟

,

R。

(

1997

)

单个tRNA的非指数动力学^PHe公司生理条件下的分子

。

程序。国家。阿卡德。科学。美国

,

94

,

7932

–

7936

。

谷歌学者

交叉参考

书目数据库

卡林

,

美国。

和

泰勒

,

H。

(

1981

)

随机过程第二课程

。

纽约

以下为：

学术出版社

。

卡林

,

美国。

和

泰勒

,

H。

(

1998

)

随机建模简介

，第3版。

纽约

以下为：

学术出版社

。

卡尔

,

答：。

(

1986

)

点过程及其统计推断

。

纽约

以下为：

德克尔

。

卡萨丁

,

R。

和

拉夫特里

,

答：。

(

1995

)

贝叶斯因子与模型不确定性

。

《美国统计杂志》。助理。

,

90

,

773

–

795

。

谷歌学者

交叉参考

书目数据库

克里切夫斯基

,

O。

和

阀盖

,

G.公司。

(

2002

)

荧光相关光谱技术及其应用

。

代表程序。物理学。

,

65

,

251

–

297

。

谷歌学者

交叉参考

书目数据库

线路接口单元

,

J.S.公司。

(

2001

)

科学计算中的蒙特卡罗策略

。

纽约

以下为：

施普林格

。

线路接口单元

,

J.S.公司。

和

萨巴蒂

,

C、。

(

2000

)

用于贝叶斯计算的广义吉布斯采样器和多重网格蒙特卡罗

。

生物特征

,

87

,

353

–

369

。

谷歌学者

交叉参考

书目数据库

卢

,

高压。

,

荀

,

L。

和

谢

,

X·S。

(

1998

)

单分子酶动力学

。

科学类

,

282

,

1877

–

1882

。

马格德

,

D。

,

埃尔森

,

E.公司。

和

韦伯

,

西。

(

1974

)

荧光相关光谱学

。

生物聚合物

,

13

,

1

–

61

。

谷歌学者

交叉参考

书目数据库

孟

,

X-L。

和

Wong（王）

,

重量小时。

(

1996

)

通过简单恒等式模拟归一化常数的比值：一种理论探索

。

统计人员。罪。

,

6

,

831

–

860

。

谷歌学者

OpenURL占位符文本

书目数据库

默纳

,

西。

(

2002

)

物理、化学和生物物理领域的十几年单分子光谱学

。

《物理学杂志》。化学。B类

,

106

,

910

–

927

。

谷歌学者

交叉参考

书目数据库

牛顿

,

文学硕士。

和

拉夫特里

,

A.E。

(

1994

)

加权似然自举的近似贝叶斯推断（附讨论）

。

J.R.统计。Soc公司。

B、，

56

,

三

–

48

。

谷歌学者

OpenURL占位符文本

书目数据库

聂

,

美国。

和

扎尔

,

R。

(

1997

)

单分子的光学检测

。

生物物理学年鉴。生物摩尔。结构。

,

26

,

567

–

596

。

谷歌学者

交叉参考

书目数据库

Pflugl公司

,

西。

,

棕色

,

F.L.公司。

和

西尔比

,

R·J。

(

1998

)

低温玻璃中单分子吸收线的变化和宽度

。

化学杂志。物理学。

,

108

,

6876

–

6883

。

谷歌学者

交叉参考

书目数据库

雷利

,

第D.页。

和

斯金纳

,

J.L公司。

(

1993

)

单分子荧光的光谱扩散：无序晶体中低频局域激发的探针

。

物理学。修订稿。

,

71

,

4257

–

4260

。

雷利

,

第D.页。

和

斯金纳

,

J.L公司。

(

1994年a

)

耦合到动态二能级系统晶格的生色团光谱：I，吸收线形状

。

化学杂志。物理学。

,

101

,

959

–

964

。

谷歌学者

交叉参考

书目数据库

雷利

,

第D.页。

和

斯金纳

,

J·L·。

(

1994年b

)

耦合到动态二能级系统晶格的生色团光谱：II，光谱扩散核

。

化学杂志。物理学。

,

101

,

965

–

973

。

谷歌学者

交叉参考

书目数据库

申特

,

G.K.公司。

,

卢

,

高压。

和

谢

,

X·S。

(

1999

)

单分子酶动力学的统计分析和理论模型

。

《物理学杂志》。化学。A类

,

103

,

10477

–

10488

。

谷歌学者

交叉参考

书目数据库

塔马拉特

,

第页。

,

马阿利

,

答：。

,

卢尼

,

B。

和

奥里特

,

M。

(

2000

)

单分子光谱学十年

。

《物理学杂志》。化学。A类

,

104

,

1

–

16

。

唐

,

J。

,

Temsamani公司

,

J。

和

阿格拉瓦尔

,

美国。

(

1993

)

自稳定反义寡核苷酸硫代磷酸酯：性质和抗HIV活性

。

核酸研究。

,

21

,

2729

–

2735

。

坦纳

,

文学硕士。

和

Wong（王）

,

重量小时。

(

1987

)

通过数据增强计算后验分布（与讨论）

。

《美国统计杂志》。助理。

,

82

,

528

–

540

。

谷歌学者

交叉参考

书目数据库

郑

,

T。

和

辛登

,

R。

(

1993

)

一级和二级DNA结构对大肠杆菌直接重复序列间缺失和重复的影响

。

遗传学

,

134

,

409

–

422

。

韦斯

,

美国。

(

2000

)

单分子荧光光谱法测量生物分子构象动力学

。

自然结构。生物。

,

7

,

724

–

729

。

谷歌学者

交叉参考

书目数据库

沃尔珀特

,

共和国。

和

艾克斯塔德

,

英国。

(

1998

)

空间统计的泊松/伽马随机场模型

。

生物特征

,

85

,

251

–

267

。

谷歌学者

交叉参考

书目数据库

谢

,

X·S。

(

2002

)

分散动力学和动力学无序的单分子方法：探索构象波动和酶动力学

。

化学杂志。物理学。

,

117

,

11024

–

11032

。

谷歌学者

OpenURL占位符文本

书目数据库

谢

,

X·S。

和

卢

,

高压。

(

1999

)

单分子酶学

。

生物学杂志。化学。

,

274

,

15967

–

15970

。

谢

,

X·S。

和

特劳特曼

,

J.K。

(

1998

)

室温下单分子的光学研究

。

A.物理版。化学。

,

49

,

441

–

480

。

谷歌学者

交叉参考

书目数据库

杨

,

H。

,

罗

,

G.公司。

,

Karnchanshamanurach公司

,

第页。

,

路易丝

,

T.-M.公司。

,

里奇

,

一、。

,

科瓦

,

美国。

,

荀

,

L。

和

谢

,

X·S。

(

2003

)

单分子电子转移法研究蛋白质构象动力学

。

科学类

,

302

,

262

–

266

。

杨

,

美国。

和

曹

,

J。

(

2001

)

单分子动力学中的双事件回声：构象起伏的特征

。

《物理学杂志》。化学。B类

,

105

,

6536

–

6549

。

谷歌学者

交叉参考

书目数据库

杨

,

美国。

和

曹

,

J。

(

2002

)

单分子动力学中记忆效应的直接测量

。

化学杂志。物理学。

,

117

,

10996

–

11009

。

谷歌学者

交叉参考

书目数据库

应

,

L。

,

华莱士

,

M。

和

克莱纳曼

,

D。

(

2001

)

发夹状DNA构象涨落的双态模型

。

化学。物理学。莱特。

,

334

,

145

–

150

。

谷歌学者

交叉参考

书目数据库

扎佐普洛斯

,

E.公司。

,

拉利

,

E.公司。

,

斯托科

,

D。

和

萨森-科尔西

,

第页。

(

1997

)

DAX-1的DNA结合和转录抑制阻断类固醇生成

。

自然

,

390

,

311

–

315

。

茨万齐格

,

R。

(

1990

)

具有动力学无序的速率过程

。

Acc.Chem.化学。物件。

,

23

,

148

–

152

。

谷歌学者

交叉参考

书目数据库

讨论中的参考文献

艾特金

,

M。

,

男孩

,

R·J。

和

查德威克

,

T·J。

(

2005

)

基于后验似然比的贝叶斯点-零假设检验

。

统计人员。计算。

，待发布。

谷歌学者

OpenURL占位符文本

书目数据库

阿斯穆森

,

美国。

(

2000

)

矩阵分析模型及其分析

。

扫描。J.统计学家。

,

27

,

193

–

226

。

谷歌学者

交叉参考

书目数据库

巴恩多夫-尼尔森

,

O.E.公司。

和

谢泼德

,

N。

(

2001

)

非高斯-奥恩斯坦-乌伦贝克模型及其在金融经济学中的一些应用（与讨论）

。

J.R.统计。Soc.B公司

,

63

,

167

–

241

。

谷歌学者

交叉参考

书目数据库

男孩

,

R·J。

,

威尔金森

,

D.J.博士。

和

柯克伍德

,

总生物量。

(

2004

)

离散观测随机动力学模型的贝叶斯推断

。

研究报告STA04.5

。

纽卡斯尔大学

,

泰恩河畔的纽卡斯尔

。

布雷莫

,

第页。

(

1981

)

点过程和队列：鞅动力学

。

纽约

以下为：

施普林格

。

布尔佐马托

,

五、。

,

比托

,

M。

,

格鲁特-科尔梅林克

,

P.J.公司。

,

科尔库洪

,

D。

和

西维洛蒂

,

L.G.公司。

(

2004

)

异聚物的单通道行为α1β甘氨酸受体：在通道开放之前检测构象变化的尝试

。

《神经科学杂志》。

,

24

,

10924

–

10940

。

卡特

,

C、。

和

科恩

,

R。

(

1996

)

条件高斯状态空间模型的马尔可夫链蒙特卡罗

。

生物特征

,

83

,

589

–

601

。

谷歌学者

交叉参考

书目数据库

Elerian公司

,

O。

,

芯片

,

美国。

和

谢泼德

,

N。

(

2001

)

离散观测非线性扩散的似然推断

。

计量经济学

,

69

,

959

–

993

。

谷歌学者

交叉参考

书目数据库

羽毛头

,

第页。

和

梅利格科茨杜

,

L。

(

2004

)

部分观测连续时间模型的精确滤波

。

J.R.统计。Soc.B公司

,

66

,

771

–

789

。

谷歌学者

交叉参考

书目数据库

羽毛头

,

第页。

和

夏洛克

,

C、。

(

2004

)

马尔可夫调制泊松过程的贝叶斯推断

。待发布。

弗雷德金

,

B.R.公司。

和

大米

,

J.A.公司。

(

1992

)

单通道膜片钳记录的贝叶斯恢复

。

生物计量学

,

48

,

427

–

448

。

格里芬

,

J·E。

和

钢材

,

M.F.J.医学博士。

(

2003

)

随机波动率的非高斯Ornstein-Uhlenbeck过程推断

。

Mimeo公司

。

华威大学

,

考文垂

。

霍克斯

,

A.G.公司。

(

2003

)单离子通道的随机建模。在

计算神经科学：一种综合方法

（编辑）。

J。

冯

)，第页。

131

–

157

。

博卡拉顿

以下为：

查普曼和霍尔——CRC

。

提花机

,

E.公司。

,

波尔森

,

N.G.公司。

和

罗西

,

体育。

(

1994

)

随机波动率模型的贝叶斯分析（附讨论）

。

J.总线。经济。统计师。

,

12

,

371

–

417

。

谷歌学者

OpenURL占位符文本

书目数据库

线路接口单元

,

J.S.公司。

和

萨巴蒂

,

C、。

(2000)

用于贝叶斯计算的广义吉布斯采样器和多重网格蒙特卡罗

。

生物特征

,

87

,

353

–

369

。

交叉参考

书目数据库

线路接口单元

,

J.S.公司。

,

Wong（王）

,

重量小时。

和

香港

,

答：。

(

1994

)

吉布斯采样器的协方差结构及其在估计量和增强方案比较中的应用

。

生物特征

,

81

,

27

–

40

。

谷歌学者

交叉参考

书目数据库

尼赫尔

,

E.公司。

、和

沙克曼

B。

(

1976

)

失神经蛙肌纤维膜记录的单通道电流

。

自然

,

260

,

799

–

802

。

牛顿

,

文学硕士。

和

拉夫特里

,

答：E。

(

1994

)

加权似然自举的近似贝叶斯推断（附讨论）

。

J.R.统计。Soc.B公司

,

56

,

三

–

48

。

谷歌学者

OpenURL占位符文本

书目数据库

帕帕斯皮利奥普洛斯

,

O。

(

2003

)

层次模型的非中心参数化和数据增强

。博士论文。

兰卡斯特大学数学与统计系

,

兰卡斯特

。

帕帕斯皮利奥普洛斯

,

O。

,

罗伯茨

,

总干事。

和

斯科尔德

,

M。

(

2003

)分层模型的非中心参数化和数据增强（带讨论）。在

贝叶斯统计7

（编辑

J·M·。

伯纳多

,

医学博士。

巴亚里

,

J.O.公司。

伯杰

,

A.第页。

Dawid公司

,

D。

赫克曼

,

上午至下午。

史密斯

和

M。

西部

)，第页。

307

–

326

。

纽约

以下为：

牛津大学出版社

。

皮特

,

M·K。

和

谢泼德

,

N。

(

1999

)

状态空间模型中Gibbs采样器的解析收敛速度和参数化问题

。

J.时间序列。分析。

,

20

,

63

–

85

。

谷歌学者

交叉参考

书目数据库

罗伯茨

,

总干事。

,

帕帕斯皮利奥普洛斯

,

O。

和

德尔拉波尔塔斯

,

第页。

(

2004

)

非高斯Ornstein–Uhlenbeck随机波动过程的贝叶斯推断

。

J.R.统计。Soc.B公司

,

66

,

369

–

393

。

谷歌学者

交叉参考

书目数据库

罗伯茨

,

总干事。

和

萨胡

,

韩国。

(

1997

)

吉布斯采样器的更新方案、相关结构、分块和参数化

。

J.R.统计。Soc.B公司

,

59

,

291

–

317

。

谷歌学者

交叉参考

书目数据库

罗伯茨

,

总干事。

和

Stramer公司

,

O。

(

2001

)

用Metropolis-Hastings算法推断部分观测非线性扩散模型

。

生物特征

,

88

,

603

–

621

。

谷歌学者

交叉参考

书目数据库

谢泼德

,

N。

和

皮特

,

M·K。

(

1997

)

非高斯测量时间序列的似然分析

。

生物特征

,

84

,

653

–

667

。

谷歌学者

交叉参考

书目数据库

威尔金森

,

D.J.博士。

,

男孩

,

R·J。

和

柯克伍德

,

总生物量。

(

2004

)

基于离散观测数据的一般随机动力学模型的贝叶斯推断

。研究报告STA04.6。

纽卡斯尔大学

,

泰恩河畔的纽卡斯尔

。

附录

附录A：数学推导

A.1、。计算 $L（左） ({t吨}_{0}, {t吨}_{1}) = P（P） {Y（Y） ({t吨}_{1}^{-}) - Y（Y） ({t吨}_{0}) = 0, γ ({t吨}_{1}) ∣ γ ({t吨}_{0})}$ 在第节中2.1

我们首先划分区间(t吨₀,t吨₁)到N个无穷小件：(t吨₀,t吨₀+小时),(t吨₀+小时,t吨₀+2小时),…,(t吨₀+(N个−1）小时,t吨₁)，其中小时= (t吨₁−t吨₀)/N个。在我第个无穷小区间是

1 - \frac{{A类}_{0}}{γ ({t吨}_{0} + 我 小时)} 小时 + o个 (小时), 我 = 0, 1, \dots, N个 - 1

这与转换矩阵一起(1.2)，表示具有离散化分辨率小时未观测到光子的概率(t吨₀,t吨₁)是

\int \prod_{我 = 0}^{N个 - 1} P（P） {(小时)}_{(我_{{t吨}_{0} + 我 小时}, 我_{{t吨}_{0} + (我 + 1) 小时})} {1 - \frac{{A类}_{0}}{γ ({t吨}_{0} + 我 小时)} 小时 + o个 (小时)} d日 γ,

（A.1）

哪里我_t吨等于1，如果γ_t吨=一和2，如果γ_t吨=b条、和P（P）(小时)_(我1,我2)表示(我₁,我₂)转移矩阵的第个项P（P）(小时). 表达式中的积分（A.1）是关于隐藏的γ过程。出租小时→0给出了概率

L（左） ({t吨}_{0}, {t吨}_{1}) = \underset{小时 \to 0}{林} [\int \prod_{我 = 0}^{N个 - 1} P（P） {(小时)}_{(我_{{t吨}_{0} + 我 小时}, 我_{{t吨}_{0} + (我 + 1) 小时)}} {1 - \frac{{A类}_{0}}{γ ({t吨}_{0} + 我 小时)} 小时 + o个 (小时)} d日 γ] 。

（A.2）

由方程式（1.2）上述乘积中的单个元素可以用简单的矩阵形式表示

P（P） {(小时)}_{(我_{0} + 我 小时, 我_{0} + (我 + 1) 小时)} {1 - \frac{{A类}_{0}}{γ ({t吨}_{0} + 我 小时)} 小时 + o个 (小时)} = {(\begin{matrix} 1 - \frac{{A类}_{0}}{一} 小时 & 0 \\ 0 & 1 - \frac{{A类}_{0}}{b条} 小时 \end{matrix}) 经验 (问 小时) + o个 (小时)}_{(我_{{t吨}_{0} + 我 小时}, 我_{0} + (我 + 1) 小时)} 。

出租

H（H） = (\begin{matrix} {A类}_{0} / 一 & 0 \\ 0 & {A类}_{0} / b条 \end{matrix}),

我们得到了方程式（A.2）作为

{[{(我 - H（H） 小时) 经验 (问 小时) + o个 (小时)}^{({t吨}_{1} - {t吨}_{0}) / 小时}]}_{{(我_{{t吨}_{0}}, 我_{{t吨}_{1}})}^{\cdot}} 。

注意到(我−H（H）小时)pt（磅）经验(Qh（质量小时）) =我+(问−H（H）)小时+o个(小时)，我们有

L（左） ({t吨}_{0}, {t吨}_{1}) = \underset{小时 \to 0}{林} {[{我 + (问 - H（H）) 小时 + o个 (小时)}^{({t吨}_{1} - {t吨}_{0}) / 小时}]}_{(我_{{t吨}_{0}}, 我_{{t吨}_{1}})} = {[经验 {(问 - H（H）) ({t吨}_{1} - {t吨}_{0})}]}_{(我_{{t吨}_{0}}, 我_{{t吨}_{1}})} 。

A.2、。表达式推导(4.10)

刘和萨巴蒂(2000)考虑了一个通用框架：比如移动一个样本x个，转换为转换后的值ν(x个)，其中转换ν是从组中提取的。他们表明，为了保持目标分布的不变性π(x个)，变换上的分布ν应该是

第页 (ν) ⑪ π {ν (x个)} | {J型}_{ν} (x个) | μ (d日 ν),

（A.3）

哪里J型(x个)是转换的雅可比矩阵μ是组上的左Haar不变测度。对于我们在这里考虑的尺度变化，左边的Haar度量μ（d）秒) =秒⁻¹d日秒雅可比式很简单秒^d日，其中秒是比例因数d日是刻度变化的尺寸。将其插入通用公式（A.3）并注意到

P（P） {秒 {x个}_{t吨} ∣ λ, \frac{{(秒 ϕ)}^{2}}{λ}} = 秒^{- (n个 + 1)} P（P） ({x个}_{t吨} ∣ λ, \frac{ϕ^{2}}{λ})

产量表达式(4.10)。

答3。定理1的证明

组移动的转换核由下式给出

\begin{array}{l} K（K） (x个, B类) = \int_{{ν 以下为： ν (x个) \in B类}} （f） (ν) 最小值 [1, \frac{π {ν (x个)} （f） ({v（v）}^{- 1})}{π (x个) （f） (ν)} | {J型}_{ν} (x个) |] μ (d日 ν) \\ + 我_{B类} (x个) (1 - \int_{Γ} （f） (ν) 最小值 [1, \frac{π {ν (x个)} （f） ({v（v）}^{- 1})}{π (x个) （f） (ν)} | {J型}_{ν} (x个) |] μ (d日 ν)) 。 \end{array}

接下来我们验证了_$X（X）$ π(x个)K（K）(x个,B类)d日x个=∫_B类π(x个)d日x个对于每个B类以便π在移动下是不变的。直接计算给出

\begin{array}{l} \int_{X（X）} π (x个) K（K） (x个, B类) d日 x个 - \int_{B类} π (x个) d日 x个 \\ \begin{array}{l} = \int_{X（X）} \int_{Γ} 我 {ν (x个) \in B类} π (x个) （f） (ν) 最小值 [1, \frac{π {ν (x个)} （f） ({v（v）}^{- 1})}{π (x个) （f） (ν)} | {J型}_{ν} (x个) |] μ (d日 ν) d日 x个 \\ - \int_{B类} \int_{T型} π (x个) （f） (ν) 最小值 [1, \frac{π {ν (x个)} （f） ({v（v）}^{- 1})}{π (x个) （f） (ν)} | {J型}_{ν} (x个) |] μ (d日 ν) d日 x个 \end{array} \end{array}

（A.4）

应用一对一映射

{\begin{array}{l} 克 = ν^{- 1}, \\ 年 = ν (x个) \end{array}

使用Fubini定理方程式（A.4）成为

\begin{array}{l} \int_{X（X）} \int_{T型} 我 {ν (x个) \in B类} π (x个) （f） (ν) 最小值 [1, \frac{π {ν (x个)} （f） ({v（v）}^{- 1})}{π (x个) （f） (ν)} | {J型}_{ν} (x个) |] μ (d日 ν) d日 x个 \\ = {¦Β}_{克 \in Γ} \int_{年 \in X（X）} 我 (年 \in B类) π {克 (年)} （f） (克^{- 1}) 最小值 [1, \frac{π (年) （f） (克)}{π {克 (年)} （f） (克^{- 1})} | \frac{\partial 年}{\partial x个} |] | \frac{\partial x个}{\partial 年} | d日 年 μ (d日 克^{- 1}) \\ \begin{matrix} = \int_{年 \in B类} \int_{克 \in Γ} 最小值 [π {克 (年)} （f） (克^{- 1}) | \frac{\partial x个}{\partial 年} |, π (年) （f） (克)] μ (d日 克^{- 1}) d日 年 \\ = \int_{B类} \int_{T型} 最小值 [π {克 (年)} （f） (克^{- 1}) | {J型}_{克} (年) |, π (年) （f） (克)] μ (d日 克^{- 1}) d日 年, \end{matrix} \end{array}

这正好是右边的最后一个术语方程式（A.4）注意到Haar测量μ（d）v（v）) =μ（d）v（v）⁻¹). 本次取消方程式（A.4）因此，这意味着

\begin{matrix} \int_{X（X）} π (x个) K（K） (x个, B类) d日 x个 = \int_{B类} π (x个) d日 x个, & 对于任何 B类 。 \end{matrix}

本文根据牛津大学出版社标准期刊出版模式的条款出版和发行(https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model)

下载所有幻灯片

月份：	总浏览次数：
2023年4月	12
2023年5月	8
2023年6月	2
2023年7月	13
2023年8月	8
2023年9月	6
2023年10月	26
2023年11月	7
2023年12月	26
2024年1月	16
2024年2月	16
2024年3月	21
2024年4月	8
2024年5月	15
2024年6月	7

文章内容

单分子实验数据的贝叶斯分析

总结

1.简介

1.1. 脱氧核糖核酸发夹

1.2. 荧光寿命实验

2.两国模型的贝叶斯分析

2.1. 可能性计算

2.2. 马尔可夫链蒙特卡罗算法

3.两态模型中的布朗扩散

3.1. 统计公式

3.2. 数据增强

3.3. 其他实验并发症

4.超越两态模型：连续扩散模型

4.1. 模型

4.2. 连续扩散模型下的似然性和数据扩充

4.3. 缩放变换更新

5.实证研究

5.1. 模拟数据集

5.2. 用两态模型分析实际数据集

5.3. 拟合连续扩散模型

5.4. 贝叶斯因子的蒙特卡罗估计

6.讨论

寇、谢和刘的论文讨论

具有潜在过程的贝叶斯推理

以前的问题

贝叶斯因子的估计

平坦先验

吉布斯采样器

贝叶斯推理、先验规范和模型敏感性

可能性计算

马尔可夫链蒙特卡罗问题

实验背景

致谢

工具书类

讨论中的参考文献

附录

附录A：数学推导

A.1、。计算L（左）(t吨0,t吨1)=P（P）{Y（Y）(t吨1−)−Y（Y）(t吨0)=0,γ(t吨1)∣γ(t吨0)}在第节中2.1

A.2、。表达式推导(4.10)

答3。定理1的证明

引文

意见

海拔高度

电子邮件警报

通过引用文章

最新的

阅读次数最多

被引用次数最多

此功能仅对订阅服务器可用

A.1、。计算 $L（左） ({t吨}_{0}, {t吨}_{1}) = P（P） {Y（Y） ({t吨}_{1}^{-}) - Y（Y） ({t吨}_{0}) = 0, γ ({t吨}_{1}) ∣ γ ({t吨}_{0})}$ 在第节中2.1