生物医学网络科学Biomedical Network Science
Example11 / 12

布尔网络的吸引子与最小钉扎控制集（穷举）

所属专题
网络可控性与稳定性
下载
Notebook · .py

6 个节点的同步更新布尔网络，先穷举 2^6 个状态找出全部吸引子，再穷举节点子集： 把子集中的节点固定为目标吸引子上的取值（钉扎），检查是否从所有初始状态都到达目标吸引子。 布尔网络的控制一般是 NP 难问题（Akutsu 等 2007），穷举只适用于很小的网络；稳定基序等方法利用结构降低搜索量。 解读：最小钉扎集通常比“全部节点”小得多，并与网络中的反馈环位置有关。

Output输出图

布尔网络的吸引子与最小钉扎控制集（穷举）：输出图 1

Output文字输出

吸引子个数 6，其中不动点：[(0, 1, 0, 0, 0, 1), (0, 1, 1, 1, 1, 0), (1, 0, 1, 1, 1, 0)]
目标不动点（A..F）： (0, 1, 1, 1, 1, 0)
最小钉扎集大小 2，共 4 种：[('A', 'C'), ('A', 'D'), ('B', 'C'), ('B', 'D')]

Code代码

# 6 个节点的同步更新布尔网络，先穷举 2^6 个状态找出全部吸引子，再穷举节点子集：
# 把子集中的节点固定为目标吸引子上的取值（钉扎），检查是否从所有初始状态都到达目标吸引子。
# 布尔网络的控制一般是 NP 难问题（Akutsu 等 2007），穷举只适用于很小的网络；稳定基序等方法利用结构降低搜索量。
# 解读：最小钉扎集通常比“全部节点”小得多，并与网络中的反馈环位置有关。
import itertools

nodes = ["A", "B", "C", "D", "E", "F"]
rules = {
    "A": lambda s: not s["B"],               # A、B 互相抑制：双稳态开关
    "B": lambda s: not s["A"],
    "C": lambda s: s["A"] or s["D"],         # C、D 互相激活：正反馈记忆环
    "D": lambda s: s["C"] and not s["F"],
    "E": lambda s: s["C"],
    "F": lambda s: s["B"] and not s["E"],
}

def step(state, pin=None):
    s = dict(zip(nodes, state))
    new = tuple(int(rules[v](s)) for v in nodes)
    if pin:
        new = tuple(pin.get(v, x) for v, x in zip(nodes, new))
    return new

def attractor_from(state, pin=None):
    seen = {}
    while state not in seen:
        seen[state] = len(seen); state = step(state, pin)
    start = seen[state]
    cyc = [s for s, k in seen.items() if k >= start]
    return frozenset(cyc)

all_states = list(itertools.product([0, 1], repeat=6))
attractors = {attractor_from(s) for s in all_states}
fixed = sorted([sorted(a)[0] for a in attractors if len(a) == 1])
print("吸引子个数 %d，其中不动点：%s" % (len(attractors), fixed))

target = max(fixed, key=lambda s: s[nodes.index("C")])     # 选一个 C=1 的不动点作为目标
print("目标不动点（A..F）：", target)
for k in range(1, 7):
    sols = []
    for S in itertools.combinations(nodes, k):
        pin = {v: target[nodes.index(v)] for v in S}
        if all(attractor_from(s, pin) == frozenset([target]) for s in all_states):
            sols.append(S)
    if sols:
        print("最小钉扎集大小 %d，共 %d 种：%s" % (k, len(sols), sols))
        break

# 作图：64 个状态的状态转移图，按所属吸引子（吸引域）着色
import networkx as nx
import matplotlib.pyplot as plt
att_list = sorted(attractors, key=lambda a: sorted(a)[0])
basin = {s: att_list.index(attractor_from(s)) for s in all_states}
STG = nx.DiGraph([(s, step(s)) for s in all_states])
fig, ax = plt.subplots(figsize=(6.4, 4.6))
pos = {}
for b, comp in enumerate(sorted(nx.weakly_connected_components(STG), key=lambda c: basin[next(iter(c))])):
    sub = nx.kamada_kawai_layout(STG.subgraph(comp).to_undirected()) if len(comp) > 1 else {next(iter(comp)): (0, 0)}
    pos.update({s: (xy[0] + 2.6 * (b % 3), xy[1] - 2.6 * (b // 3)) for s, xy in sub.items()})
nx.draw_networkx_edges(STG, pos, ax=ax, edge_color="0.6", alpha=0.5, arrowsize=6, node_size=40)
cmap = plt.get_cmap("tab10")
for b, a in enumerate(att_list):
    members = [s for s in all_states if basin[s] == b]
    lab = "attractor %d (%s, %d states)" % (b, "fixed point" if len(a) == 1 else "cycle %d" % len(a), len(members))
    nx.draw_networkx_nodes(STG, pos, nodelist=members, ax=ax, node_size=[120 if s in a else 40 for s in members],
                           node_color=[cmap(b)], label=lab)
ax.legend(fontsize=7, loc="upper center", bbox_to_anchor=(0.5, -0.02), ncol=2); ax.axis("off")
ax.set_title("State transition graph (large dots = attractor states)", fontsize=10)
plt.tight_layout(); plt.show()

Run在浏览器中运行

首次加载 Python 环境约需 20–40 秒；可以修改参数后重新运行。

布尔网络的吸引子与最小钉扎控制集（穷举）