Python莫兰生死抑制放大进化图

🎯要点

🎯种群离散时间莫兰生死动态图解 | 🎯良好混合种群的固定概率 | 🎯数值求解生成埃尔多斯-雷尼图 | 🎯计算马尔可夫链的转移矩阵概率 | 🎯出生死亡动态和死亡出生动态概率无向随机图和有向随机图，转移矩阵概率计算

📜进化图用例

📜Python种群邻接矩阵彗星风筝进化图算法

📜Python和C++骨髓细胞进化解析数学模型

🍪语言内容分比

🍇Python转移矩阵

令 $E$ 定义构成时间序列数据的 $k$ 个唯一事件集。例如，时间序列可能由以下三个基本且唯一的事件组成，这些事件表示在离散时间步长上绘制数据时观察到的路径轨迹类型：向下、持平和向上。令$S$ 定义长度为 $n$（表示离散时间步长）的序列，该序列由$E$ 中定义的事件组成，表示部分或全部数据。例如，序列 [向上、向下、向上、持平、向上] 表示五个时间步长的数据。

现在可以定义一个维度为 $k^2$的马尔可夫转移矩阵 $M$，使得每个元素 $M(i,j)$ 描述在给定时间序列中从时间步骤 $t$ 中的事件 $E(i)$转移到时间步骤 $t+1$ 中的事件 $E(j)$ 的概率。换句话说，$M(i, j) $表示在连续时间步骤中两个事件之间转移的条件概率。从图论意义上讲，如果时间序列数据中 $E(i)$ 后面是 $E(j)$，则事件$ E(i)$ 和 $E(j)$ 可以被认为是由有向边 $E(i)\to E(j)$ 连接的节点，那么马尔可夫转移矩阵 $M$​​ 本质上表示图中节点所描绘事件的邻接矩阵（或共现矩阵）的规范化版本。

假设我们有以下涵盖 11 个连续时间步长的原始时间序列数据：[1, 2, -2, -1, 0, 0, 2, 2, 1, 2, 3]。使用上面描述的路径轨迹的简化视图，我们可以将数据转换为以下 10 个事件序列，这些事件描述相邻时间步长之间的转换：[上、下、上、上、平、上、平、下、上、上]。

我们现在可以构建以下邻接矩阵来捕获事件序列中同时出现的模式：
$$A=\left(\begin{array}{lll}2& 2& 1\\ 1& 0& 1\\ 2& 0& 0\end{array}\right)$$
元素 $A(i,j)$ 表示事件序列中某个时间步 $t$ 的事件 $i$ 后面跟着时间步 $t+1$ 的事件 $j$ 的次数； $i$ 和 $j$ 分别是行索引和列索引。请注意，行表示从上到下、从上到下的顺序的事件，而列从左到右表示相同的事件。例如，$A$ 的左上角元素表示在给定的事件序列中，上事件后紧跟着另一个上事件两次。 $A$ 的中右元素表示在事件序列中，平事件之后紧接着下事件。等等。

我们可以按行或按列标准化矩阵 $A$ 以生成转换矩阵。如果我们使用基于行的归一化，则元素 $M(i,j)$ 将描述给定事件 $E(i)$ 在时间步 $t+1$ 中看到事件 $E(j)$ 的概率时间步$t$。因此，每行中的概率之和应为 1 。在我们的示例中，行归一化矩阵如下所示：
$$M_{\text{rnorm }}=\left(\begin{array}{ccc}2/5& 2/5& 1/5\\ 1/2& 0& 1/2\\ 2/2& 0& 0\end{array}\right)=\left(\begin{array}{ccc}0.4& 0.4& 0.2\\ 0.5& 0& 0.5\\ 1& 0& 0\end{array}\right)$$
类似地，如果我们要使用基于列的归一化，则元素 $M(i,j)$ 将描述给定时间步 $t$ 中的事件 $E(j)$ 的情况下，在时间步 $t-1$ 中发生事件 $E(i)$ 的概率。现在每列中的概率之和应为 1。在我们的示例中，列归一化矩阵如下所示：
$$M_{\text{cnorm }}=\left(\begin{array}{ccc}2/5& 2/2& 1/2\\ 1/5& 0& 1/2\\ 2/5& 0& 0\end{array}\right)=\left(\begin{array}{ccc}0.4& 1& 0.5\\ 0.2& 0& 0.5\\ 0.4& 0& 0\end{array}\right)$$
请注意，行归一化（名义上向前看时间）的条件概率可能与列归一化（向后看时间）的条件概率不同。

代码计算

import pandas as pddef get_transition_tuples(ls):return [(ls[i-1], ls[i]) for i in range(1, len(ls))]def get_transition_event(tup):transition_event = 'flat'if tup[0] < tup[1]:transition_event = 'up'if tup[0] > tup[1]:transition_event = 'down'return transition_eventls_raw_time_series = [1, 2, -2, -1, 0, 0, 2, 2, 1, 2, 3]
ls_transitions = get_transition_tuples(ls_raw_time_series)
ls_events = [get_transition_event(tup) for tup in ls_transitions]
ls_event_transitions = get_transition_tuples(ls_events)
ls_index = ['up', 'flat', 'down']
df = pd.DataFrame(0, index=ls_index, columns=ls_index)for i, j in ls_event_transitions:df[j][i] += 1  df_rnorm = df.div(df.sum(axis=1), axis=0).fillna(0.00)df_cnorm = df.div(df.sum(axis=0), axis=1).fillna(0.00)

这应该产生以下转换矩阵：

>>> df  up    flat  down
up    2     2     1
flat  1     0     1
down  2     0     0
>>> df_rnorm  up    flat  down
up    0.4   0.4   0.2
flat  0.5   0.0   0.5
down  1.0   0.0   0.0
>>> df_cnorm  up    flat  down
up    0.4   1.0   0.5
flat  0.2   0.0   0.5
down  0.4   0.0   0.0

可视化转移状态

def get_df_edgelist(df, ls_index):edgelist = []for i in ls_index:for j in ls_index:edgelist.append([i, j, df[j][i]])return pd.DataFrame(edgelist, columns=['src', 'dst', 'weight'])def edgelist_to_digraph(df_edgelist):g = Digraph(format='jpeg')g.attr(rankdir='LR', size='30')g.attr('node', shape='circle')nodelist = []for _, row in df_edgelist.iterrows():node1, node2, weight = [str(item) for item in row]if node1 not in nodelist:g.node(node1, **{'width': '1', 'height': '1'})nodelist.append(node1)if node2 not in nodelist:g.node(node2, **{'width': '1', 'height': '1'})nodelist.append(node2)g.edge(node1, node2, label=weight)return gdef render_graph(fname, df, ls_index):df_edgelist = get_df_edgelist(df, ls_index)g = edgelist_to_digraph(df_edgelist)g.render(fname, view=True)

render_graph('adjmat', df, ls_index)
render_graph('transmat_rnorm', df_rnorm, ls_index)
render_graph('transmat_cnorm', df_cnorm, ls_index)

👉参阅、更新：计算思维 | 亚图跨际