文献复现

发布于 2023-07-11  320 次阅读


Please refresh the page if equations are not rendered correctly.
---------------------------------------------------------------

3.2.1概率转移矩阵的规范形式

寻求一个概率转移矩阵来再现实验数据的两个统计量 RMSD \sigma_{\xi}^{(w)}和相关长度 \zeta_{\xi}^{(w)}(\xi\equiv y\mathrm{~or~}z) (经纱的例子)。这两个数据是通过由平稳马尔科夫算子A产生的概率分布的渐进行为来实现拟合的。

A多次迭代作用(重复应用P_{i}^{(n+1)}=AP_{i}^{(n)}),得到的分布向量逐渐地接近与A的最大特征值相关的A的特征向量(任何初始分布向量都可以写成A的特征向量p_{k}的线性组合,A,k=1,\ldots,2m+1P(n)=A^{n}P(0)=\sum_{k=1}^{t}C_{k}\mu_{k}^{n}p_{k}。当A乘以初始向量,线性组合中第K个特征向量的系数乘以相应的特征值\mu_{k},如果重复n次,则每个系数相乘的因子为\mu_{k}^{n}。所以对于较大的n\mu_{k}最大的项占主导地位。)设\mu^{>}表示最大特征值,p^{>}表示相应的特征向量,因为A的元素都是正的,列和为1,\mu^{>}=1。

证:为什么达到稳态时,\mu^{>}=1。

由于A的每一列和都为1,所以考虑用A-E ,显然有:

|A-E|=0

A-E相当于A的每一列和都减去1,因此A-E的每一列和都为0。

|A-E|=\prod_{k=1}^{2m+1}(\mu_{k}-1)

(矩阵特征值之积等于矩阵行列式的值)

于是有

\prod_{k=1}^{2m+1}(\mu_{k-1})=0

所以一定存在特征值为1的情况。且A的列和为1,所以

|\mu_{k}|\leq1

A的最大特征值为1,得证。

Markov chain实现的 RMSD 可以用p^{>}表示:

\tilde{\sigma}_{\xi}^{(w)}= \sqrt{ \langle \delta \tilde{\xi}^{2} \rangle } \newline =\left[ \frac{\sum_{k}{(ka)}^{2}\Pr \lbrace \delta \tilde{\xi}^{(n,j,w)}=ka \rbrace }{\sum_{k}\Pr \lbrace \delta \tilde{\xi}^{(n,j,w)}=ka \rbrace } \right]^ {1/2}

\langle \rangle表示期望值,\Pr \lbrace {\delta\tilde{\xi}^{(n,j,w)}=ka} \rbrace\delta\tilde{\xi}^{(n,j,w)} 在网格点 (n,j) 取离散值 ka 的概率。其中分母\sum_{k}\mathrm{Pr}{\delta\tilde{\xi}^{(n,j,w)}=ka}表示在网格点(n,j)取离散值ka的概率的全部加和,概率和为1。

在马尔科夫过程的稳态状态中,\mathrm{Pr}{\delta \tilde{\xi}^{(n,j,w)} = ka}=p_{k}^{>}p^{>}中的第k个元素)有

\tilde{\sigma}_\xi^{(w)}=a\sqrt{\sum_kk^2p_k^>}

在连续的数据点\delta\tilde{\xi}^{(n,j,w)}的值之间实现的相关性为

\tilde{\mathcal{C}}_\xi^{(n,n+1,j,w)}=\langle\delta\tilde{\xi}^{(n,j,w)}\delta\tilde{\xi}^{(n+1,j,w)}\rangle/(\tilde{\sigma}_\xi^{(w)})^2=\left(\frac{1}{\tilde{\sigma}_\xi^{(w)}}\right)^2\frac{\sum_{k,k}kak^{\prime}a\Pr \lbrace(\delta\tilde{\xi}^{(n,j,w)}=ka)\&(\delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a) \rbrace }{\sum_{k,k}\Pr \lbrace (\delta\tilde{\xi}^{(n,j,w)}=ka)\&(\delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a) \rbrace}

\Pr \lbrace (\delta\tilde{\xi}^{(n,j,w)}=ka)\&(\delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a) \rbrace\delta\tilde{\xi}^{(n,j,w)}在沿丝束的连续网格点取离散值kak^{\prime}a的联合概率。

在稳定状态下

\Pr \lbrace (\delta\tilde{\xi}^{(n,j,w)}=ka)\&(\delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a) \rbrace =\Pr \lbrace \delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a|\delta\tilde{\xi}^{(n,j,w)}=ka \rbrace \Pr \lbrace \delta\tilde{\xi}^{(n,j,w)}=ka \rbrace =A_{k^{\prime},k}p_k^{>}

所以

\sum_{k,k^{\prime}}\Pr \lbrace (\delta\tilde{\xi}^{(n,j,w)}=ka)\&(\delta\tilde{\xi}^{(n+1,j,w)}=k^{\prime}a) \rbrace =\mu^{>}\sum_{k}p_{k}^{>}=1

(因为Ap_{k}=\mu_{k}p_{k},在矩阵乘法中,矩阵乘以特征向量等于对应的特征值乘以特征向量,\mu^{>}=1,\sum_{k}p_{k}^{>}=1p_{k}^{>}一组归一化的概率也为1)。

\tilde{C}_\xi^{(n,n+1,j,w)}=\frac{\langle\delta\tilde{\xi}^{(n,j,w)}\delta\tilde{\xi}^{(n+1,j,w)}\rangle}{(\tilde{\sigma}_\xi^{(w)})^2}=\left(\frac{a}{\tilde{\sigma}_\xi^{(w)}}\right)^2\sum_{k,k}kk^{\prime}p_k^>A_{k^{\prime},k}

根据公式的假设,其独立于网格数n或纤维束j。所以现在考虑(2m+1)\times(2m+1)三对角概率转移矩阵。

A_{0} =\begin{vmatrix} \alpha& \gamma & 0& \cdot&&&&&& \cdot& 0 \\ {\beta}'& \alpha& \gamma& 0& \cdot& &&&&& \cdot \\ 0& \beta& \alpha & \cdot& \cdot& \cdot& &&&& \\ \cdot& 0& \cdot& \cdot& \gamma& 0& \cdot& &&&\\ & \cdot& \cdot& \beta& \alpha& \delta& 0& \cdot& &&\\ && \cdot& 0& \beta& \alpha& \beta& 0& \cdot& &\\ &&& \cdot& 0& \delta& \alpha& \beta& \cdot& \cdot& \\ &&&& \cdot& 0& \gamma& \cdot& \cdot& 0& \cdot \\ &&&&& \cdot& \cdot& \cdot& \alpha& \beta&0 \\ \cdot&&&&&&& 0& \gamma& \alpha&{\beta}' \\ 0& \cdot& &&&&& \cdot& 0& \gamma&\alpha \end{vmatrix}

其中所有元素取[0,1]中的值,所有列相加为一。即 α+β+γ=1α+2δ=1α+\beta^{\prime}=1。且A_{0}具有反转对称性(A_{k,k^{\prime}}=A_{-k,-k^{\prime}},0\leq k\leq m,0\leq k^{\prime}\leq m)

可以把偏差\delta\tilde{\xi}^{(n,j,w)}离散到\lbrace -ma,-(m-1)a,..,0,...,(m-1)a,ma \rbrace,\mathrm{~with~}ma=3\sigma_\xi^{(w)}看做是一个生死链过程(birth-death chain),即

如果A_{0}被用作生成器,则\delta \tilde{\xi}^{(n,j,w)}将在零附近均匀分布。所产生的偏差 \tilde{\sigma}_{\xi}(省略了丝束两种类型符号)可以通过改变\beta\gamma的相对大小来控制。当\beta大于\gamma时,在birth-death chain的每一步都倾向于向\delta\tilde{\xi}=0靠拢,所以当\beta - \gamma从0开始增加时,\tilde{\sigma}_{\xi}减小。

在极限β=γ时,β=γ=δ,birth-death chain变为

如果不考虑边缘离散值m和-m的概率,可以近似看做每个网格点离散到每一个相应值的个数都是不变的且都是相等的,所以\tilde{\sigma}_{i}可以取以下近似值

\sigma_{max}=\sqrt{\frac{\sum_{k=-m}^m{(ma)}^2}{2m+1}}=\sqrt{\frac{2\sum_{k=0}^m{m}^2}{2m+1}}a=\sqrt{\frac{2\cdot\frac{1}{6}(m+1)(2m+1)m}{2m+1}}a=\sqrt{\frac{m(m+1)}3}a

1^{2}+2^{2}+3^{2}+\cdots+m^{2}=\frac{m(m+1)(m+2)}{6}

\tilde{\sigma}_{i}近似等于\tilde{\sigma}_{max}的前提是忽略了边缘离散值m和-m的概率,所以\tilde{\sigma}_{max}\tilde{\sigma}_{i}的极限值之间的差异就取决于边缘离散值的概率与中间离散值概率值比β/\beta^{\prime},图3a显示了\tilde{\sigma}_{\xi}/\sigma_{\max}的变化,如果 α变化,β也变化以保持β/\beta^{\prime}也就是\beta/(1-\alpha)不变那么\tilde{\sigma}_{\xi}/\sigma_{\max}也保持不变。\tilde{\sigma}_{\xi}/\sigma_{\max}随m的变化不显著,可以通过选择适当的ma,找到β来达到任何目标的RMSD。

无论 αβ-γ的值是多少,相关参数\tilde{C}_{\xi}^{(n,n+1)}在数值上保持接近1(忽略对丝束类型的引用)。在这里我的理解是参考文章前面的假设(a3)同一丝束上不同点处的分量 \delta\tilde{\xi} 是相关的,相关程度随着点的分离而减弱。相邻两个网格点的相关程度很大,所以相关参数在数值上保持接近1.如果\alpha\rightarrow1\tilde{C}_\xi^{(n,n+1)}\to1(也就是在birth-death chain中,不同网格点的离散值保持不变),马尔科夫序列中\delta\tilde{\xi}的连续值是完全相关的。随着α的减少,\tilde{C}_{\xi}^{(n,n+1)}的值略微减小,但是 \delta\tilde{\xi} 连续值的间隔永远不超过a,因此\tilde{C}_{\xi}^{(n,n+1)}减小会受到限制。A_{0}的三对角形式可以生成任何的RMSD分布,但是不能独立控制相关参数。

因此,考虑概率转移矩阵

A_N=A_0^N

其表示三对角转移矩阵A_{0}N个迭代作用,首先要注意的是

A_Np^>=A_0p^>

因为在稳定状态下

A_{0}P^{>}=A_{0}^{N}P^{>}=A_{N}P^{>}

即在稳定状态下,A_{0}A_{N}产生相同的 \delta\tilde{\xi} 分布。因此,如果选择A_{0}来实现\tilde{\sigma}_{\xi}=\sigma_{\xi}^{(w)} ,那么A_{N}也将实现\tilde{\sigma}_{\xi}=\sigma_{\xi}^{(w)} 。但是如果不是在稳定状态下,相关参数\tilde{C}_{\xi}^{(n,n+1)}将随着N的增加而单调减小。因为马尔科夫链中初始状态分布必须随着每次过渡到新状态而逐渐消失,即P_{i}^{(n)}\rightarrow A_{0}P_{i}^{(n)}。即使α非常接近1,如果N选择足够大, \tilde{C}_{\xi}^{(n,n+1)}也可以减小到零。

3.2.2校准马尔科夫链

首先为离散变量\delta\xi创建网格,要求ma满足

\sigma_\xi=\frac{ma}3

即允许在实验确定的标准偏差上波动,波动不能超过3个偏差 ,从而使得离散空间\lbrace -ma,-(m-1)a,\ldots,0,\ldots,(m-1)a,ma \rbrace可以代表\delta\xi的值。

选择m=10,与最可能的数值范围相比,可以确保\delta\xi在所允许的离散值之间的间距是相当小的,那么a的值就可由上式得出。

如果选择α=0.9,则条件

\tilde{\sigma}_{\xi}=\sigma_{\xi}

减小到

f=\frac{m/3}{\sqrt{\frac13m(1+m)}}\approx0.55\quad\text{(for }m=10)

f(\beta,\alpha)是图a中的函数,这与\beta+\gamma=1-\alpha一起决定了\beta的值。在设定了\beta的情况下,通过要求生成的相关参数\tilde{C}_{\xi}^{(n,n+1)}等于实验确定的统计量C_{\xi}来确定N(如图b)。

用论文前面2.4中的数据,对每个丝束属重复校准。

3.2.3操作虚拟样本生成器

虚拟样本生成器通过将质心坐标 (\xi,z) 分解为非随机的、周期性变化(系统性变化),再加上随机性的偏差。

(y,z)^{(n,j,w)}=\langle(y,z)^{(n,j,w)}\rangle+(\delta y,\delta z)^{(n,j,w)}

对坐标的系统性变化是通过应用以下方程产生的。

\langle(y,z)^{(nj,w)}\rangle=\boldsymbol{u}_\phi^{(w)}+r\frac{\lambda_x}4\hat{\boldsymbol{x}}+s\frac{\lambda_y}4\hat{\boldsymbol{y}}

通过Monte Carlo生成每根丝束的随机变化。生成器生成的随机坐标的统计数据对于任何丝束上的所有统计点都应该是统一的。在纤维的开始处不应该存在马尔科夫过程的结果逐渐向稳态分布p^{>}过渡,因为Monte Carlo是基于马尔科夫链平稳分布的对应样本集进行采样的。可以通过从分布p^{>}中选择丝束上第一个点的随机坐标值来保证样本集均是平稳分布的。当在许多副本上采样时,通过应用概率转移矩阵为丝束上的后续点生成的坐标值集合将具有相同的分布。

Monte Carlo方法的实现涉及以下从概率分布向量导出的分段常数函数,

3.3消除网格效果

马尔科夫链算法以与目标RMSD和相关长度一致的方式成功地再现了沿纤维束轴变化的随机坐标偏差。

图中显示了随机偏差(中心纬纱丝束的)的典型实例。大振幅、长波长变化发生在沿丝束的距离上,该距离与目标相关长度\zeta_{x}^{(cw)}相当。然而,该马尔科夫链算法也会产生小幅度、短波长的振荡,这些振荡在原始生成的曲线上表现为尖锐的随机尖峰(也就是图中的细线)。尖峰的振幅可以通过细化随机变量的网格来减小(增加离散值的m,但是要同时保持ma=3\sigma_{x}),并且尖峰的波长取决于沿丝束的网格密度。因此尖峰是依赖于网格的。

但是这些尖峰是不现实的,因为连接的丝束不可能在其位置上表现出如此强烈的变化。尖峰是由于假设偏差被视为马尔科夫过程计算产生的。在图中的例子,网格点k处的值\delta x_{k}仅取决于网格点k-1处的\delta x_{k-1},而不取决于网格点k-2处的\delta x_{k-2}。马尔科夫链生成器没有模仿连续纤维保持平滑倾向的约束。

可以通过对生成的偏差序列应用局部平均来平滑尖峰。为了最大限度地减少其对生成偏差的RMSD的影响,使用新制定的算法对偏差的平方进行平均。

\chi_{i}^{(s)}=S\sqrt{\frac{1}{2m+1}\sum_{j=-m}^m\chi_{i+j}^2}

上标s代表平滑的量, \chi表示偏差,等式中的S是符号因子,因为偏差有正负,而取平方根的结果是非负的,所以取S来表示偏差的正负。该等式仅适用于平均区间内所有\chi_{i}具有相同的符号或者为零的情况,非零项的符号就是S的符号。当平均区间内穿过\chi=0时,该等式不适用。

对于 \chi_{i}中包含符号变化的区间,使用公式常规平均规则。

\chi_{i}^{(s)}=\frac{1}{2m+1}\sum_{j=-m}^m\chi_{i+j}

等式\chi_{i}^{(s)}=S\sqrt{\frac{1}{2m+1}\sum_{j=-m}^m\chi_{i+j}^2}在偶尔出现平均区间\chi=0穿过零的情况下也能最大程度地减少平滑对生成偏差的RMSD的影响。在马尔科夫链算法产生的偏差中,穿过零的情况相对不多,因为与偏差的相关长度相比,网格是密集的,而且在零附近的元素对RMSD的影响很小。

为了避免影响丝束的相关性,在小于偏差的基准相关长度的标距长度进行平均。细平滑曲线是对偏差序列应用常规平均规则平滑的结果,粗平滑曲线是应用新制定的算法对偏差的平方进行平滑的结果。对于每个平滑过程,使用±4个网格间隔的加权信息递归平滑操作,从而合并了在任意网格点i到网格点i±4的加权信息。

粗曲线所示:物理上显著的大振幅、长波长变化保持不变,而非物理的尖峰则减少到振幅非常小(几毫米)的波纹,对复合材料的预测机械性能的影响可以忽略不计。

届ける言葉を今は育ててる
最后更新于 2023-07-11