简介：基于Jupyter Notebook的生成对抗网络的Hi-C数据增强方法，面向毕业设计、课程设计及项目开发场景，适合具备一定Python与深度学习基础的读者参考。包体共92个文件，以ipynb交互式笔记本、py脚本、模型权重（pth）、收敛日志（xls）及可视化图片（png）为主，压缩包约173.27MB；内含标准GAN、CGAN、WGAN、WGAN-GP等多种训练实现，并划分出utils、xxydataset、models等模块，分别封装数据处理、数据集与生成器/判别器定义。项目已通过严格测试，附带收敛曲线记录、模型距离测试、最优模型筛选及图像生成与展示脚本，目录结构清晰，便于按需复现与二次开发。目前已有97人学习下载，可作为Hi-C数据增强相关课题的完整参考实现。
1. Hi-C 数据增强为什么需要对抗生成网络
把“生成对抗网络”塞进 Hi-C 数据增强任务，十个人里有九个会先去找现成的 DCGAN 代码换数据跑，结果出来的矩阵要么一坨噪声，要么就是把上采样插值锐化了一遍。在 Jupyter Notebook 里真正决定成败的从来不是网络结构，而是 Hi-C 接触矩阵本身的三个特征：样本少、稀疏度高、沿对角线的距离衰减曲线极强，它们决定了训练对怎么构造、损失函数怎么写。下面这套方案是课程设计最常见的落地形态：从 .cool 文件读矩阵，切块构造低分辨率到高分辨率的训练对，用 PyTorch 在 Notebook 里完成对抗生成网络训练，最后给出可验证指标和排错清单。适合正在写期末大作业的学生，也适合需要扩增三维基因组数据集的工程岗。
2. Hi-C 接触矩阵的预处理：从 .cool 到可训练的张量
装 Jupyter 不是重点，按安装教程用 conda 建一个 py39 环境再 pip install notebook 是最省事的装法；进入主题之后，第一件事是明确 Hi-C 数据增强的输入形态。大多数公开数据集的接触矩阵来自 Hi-C 实验，存储格式要么是 .hic，要么是 cooler 生态里的 .cool。你不需要从原始测序 reads 自己跑比对和分 bin，常见做法是先用 hic2cool 或 hicConvertFormat 把 .hic 转成 .cool，再在 Notebook 里只负责矩阵层面的加工。这一步做完，后面 GAN 跑得快不快、生成结果像不像，80% 都在这一步定下来。
2.1 选任务：低分辨率重建还是合成样本生成
数据增强在 Hi-C 场景下主要有两条路线，很多人一开始会混着做，导致训练目标和评价指标对不上。
| 任务类型 | 输入 | 输出 | 验证难度 | 适用场景 |
|---|---|---|---|---|
| 低分辨率重建 | 下采样后的低分辨率矩阵 | 原始分辨率矩阵 | 低，有 ground truth，可直接算 SSIM/PCC | 修复稀疏测序样本、补全降采样数据 |
| 合成样本生成 | 随机噪声向量 | 与真实矩阵相似的完整矩阵 | 高，只能靠分布距离和结构特征间接验证 | 扩充下游模型训练集、异常检测训练 |
课程设计和期末大作业我一般建议选低分辨率重建。原因有两个：一是评价指标客观，老师可以拿着真实矩阵对答案；二是这个任务能复用公开的染色质构象捕获数据，不需要额外造标签。合成样本生成当扩展练习做，不要当主任务。
2.1.1 下采样方式必须和 Hi-C 的分辨率语义一致
天然 Hi-C 分辨率提升，是把多个小 bin 合并成大 bin，接触计数求和或取平均。用这个语义做训练对，模型学到的才是“分辨率升级”，而不是单纯把图放大。代码里我用 mean pooling 模拟 bin 合并，再通过线性插值把尺寸还原回原图大小，作为生成器的条件输入。
2.2 读取矩阵并切块：最小可运行代码
先用 cooler 把某个染色体区间读成二维矩阵。如果手头是远程集群上的 .hic 文件，可以先用命令行导出文本矩阵，再在 Notebook 里用 numpy 装载，这里给的是最常用的一条路。
import numpy as np import cooler from scipy.ndimage import zoom # 读取 50kb 分辨率的 cool 文件，balance=True 拿 ICE 归一化后的矩阵 uri = "HiC_50kb.cool" clr = cooler.Cooler(uri) region = "chr7:10000000-20000000" # 10 Mb 窗口，50kb 下得到 200 x 200 raw = clr.matrix(balance=True).fetch(region) mat = np.nan_to_num(raw) # 未覆盖区域是 NaN，直接填 0 print("matrix shape:", mat.shape) # (200, 200)参数说明：balance=True会返回做过归一化的矩阵，空白区和低覆盖 bin 会是 NaN，nan_to_num统一填成 0，防止后面切块时把 NaN 带进 GAN 输入。分辨率的选择直接影响矩阵大小，10 Mb 窗口配 50kb bin 是 200x200，配 25kb 就是 400x400，后者对显存很不友好，起步阶段先用 50kb。
接下来做切块。整张矩阵直接送进网络在显存上不现实，标准做法是按固定窗口重叠切块：
def to_patches(mat, size=64, stride=32, q=0.99): h, w = mat.shape if h < size: mat = np.pad(mat, ((0, size - h), (0, size - w))) patches = [] coords = [] for i in range(0, h - size + 1, stride): for j in range(0, w - size + 1, stride): block = mat[i:i + size, j:j + size] vmax = np.quantile(block, q) # 用 99 分位砍掉极端接触 block = np.clip(block / vmax, 0.0, 1.0) patches.append(block) coords.append((i, j)) return np.stack(patches)[..., None], np.array(coords)逻辑说明：每个 patch 单独做分位数归一化，把数值压到 0-1，避免基因组上不同区域接触强度跨度太大导致 GAN 训练震荡。q=0.99的意思是只把超过 99 分位的极端接触截断，那些通常是核仁或中心粒造成的伪峰。stride 取 size 的一半形成 50% 重叠，数据量翻倍，模型能看到更多局部上下文。
2.3 构造低分辨率训练对的函数
有了高分辨率的 patch，低分辨率输入按 bin 合并语义下采样，再插值回原尺寸：
def make_pair(high, factor=2): high = high.squeeze().astype(np.float32) n = high.shape[0] // factor sub = high[:n * factor, :n * factor] low_raw = sub.reshape(n, factor, n, factor).mean(axis=(1, 3)) # bin 合并 low = zoom(low_raw, factor, order=1) # 插值回原尺寸 return low[..., None], high[..., None]参数说明：factor=2表示输入的低分辨率矩阵是目标的一半分辨率；把这个值提高到 4 可以让任务更难，但生成器需要更大感受野。注意这里不用cv2.resize做双线性插值，因为scipy.ndimage.zoom对二维数组更直接，order=1是线性插值，生成器要学的是线性插值结果与真实高分辨率之间的残差，这比从纯几何上采样重建容易得多。
代码里用mean而不是sum做 bin 合并，是为了与zoom的量纲保持一致，也避免大 bin 里接触计数天然更大而把输入的整体幅值放大。最后把训练集按染色体拆分，绝对不要将同一条染色体的相邻区域同时放进 train 和 validation，否则评估结果会因为空间自相关虚高。切完的大数组直接转成torch.Tensor并放进TensorDataset，比每个 epoch 现读文件快得多。
3. 用 PyTorch 在 Jupyter Notebook 里搭一个 Hi-C GAN
模型结构不需要套用 ImageNet 上那套复杂的设计，但也不能直接抄 DCGAN。Hi-C 矩阵和自然图像有个根本差异：真实的接触信号集中在对角线区域，远离对角线的区域几乎全是接近 0 的低值。DCGAN 的转置卷积堆叠很容易把注意力放在高频纹理上，生成一堆角落噪声，而结构上该保留的大尺度 compartment 反而丢了。
3.1 生成器：U-Net 式残差结构
输入输出尺寸相同，所以生成器采用编码器-解码器加残差的组合。核心思路是让网络学习“低分辨率输入到高分辨率输出之间的误差图”，而不是凭空生成矩阵。
import torch import torch.nn as nn class ResBlock(nn.Module): def __init__(self, ch): super().__init__() self.net = nn.Sequential( nn.InstanceNorm2d(ch), nn.LeakyReLU(0.2, inplace=True), nn.Conv2d(ch, ch, 3, 1, 1), nn.InstanceNorm2d(ch), nn.LeakyReLU(0.2, inplace=True), nn.Conv2d(ch, ch, 3, 1, 1), ) def forward(self, x): return self.net(x) + x class Generator(nn.Module): def __init__(self, base=32): super().__init__() self.head = nn.Conv2d(1, base, 3, 1, 1) self.enc = nn.Sequential( ResBlock(base), nn.Conv2d(base, base * 2, 4, 2, 1), ResBlock(base * 2), ) self.dec = nn.Sequential( nn.ConvTranspose2d(base * 2, base, 4, 2, 1), ResBlock(base), nn.Conv2d(base, 1, 3, 1, 1), ) def forward(self, x): return self.dec(self.enc(self.head(x))) + x这里的逻辑是：ResBlock内部用 InstanceNorm 而不是 BatchNorm，是因为每个 batch 里的 patch 来自基因组上不同区域，接触强度分布差异很大，BN 会把这种差异当噪声抹掉。外层self.dec(...) + x是一个全局残差连接，网络只需要输出残余的细节，训练初期 Loss 就能在较低水平起步，不至于从空白图开始学。
3.1.1 为什么不加深编码器
64x64 的 patch 用两层下采样到 16x16 已经足够捕捉 compartment 级别的结构。再往下采样到 8x8 并不会带来明显收益，反而会拖慢训练。真正影响生成质量的是训练对构造时的factor，它决定生成器要补的细节量级。
3.2 判别器：PatchGAN 而不是全局标量
判别器输出一个局部概率图，让每个感受野都独立判断它看到的矩阵片断像不像真实 Hi-C patch，比全局判别器更能捕捉对角线附近的局部结构。
class Discriminator(nn.Module): def __init__(self, base=32): super().__init__() self.net = nn.Sequential( nn.Conv2d(1, base, 4, 2, 1), nn.LeakyReLU(0.2, inplace=True), nn.Conv2d(base, base * 2, 4, 2, 1), nn.InstanceNorm2d(base * 2), nn.LeakyReLU(0.2, inplace=True), nn.Conv2d(base * 2, 1, 4, 1, 1), ) def forward(self, x): return self.net(x)输出 shape 是(batch, 1, 8, 8)，对应原图上 8x8 个重叠感受野块。损失计算时，每个位置都算一次 BCE，再取平均。这样即使某个区域的少数 patch 造假得特别好，也只是拉低该位置的分数，不至于让整个样本通过。
3.3 训练循环和 Hi-C 专用的距离衰减损失
对抗损失用 BCE，同时在生成器上加两个监督项：L1 保持整体强度，再加一个按对角线层级计算的损失，让每个线性距离的接触强度都被同等地监督到。
def decay_loss(pred, target): diff = (pred - target).abs() n = diff.shape[-1] total = 0.0 for k in range(n): diag = diff.diagonal(offset=k, dim1=-2, dim2=-1) total = total + diag.mean() return total / n逻辑说明：对角线上第k个元素对应线性距离k个 bin 的接触对。把每条对角线的平均绝对误差单独算再平均，等价于对每个距离层级给出相等的权重。如果只用全局 L1，模型会优先拟合接触值大的短距离区域，长距离的远端相互作用（compartment 信号所在的位置）被直接忽略。
def train_one_epoch(G, D, loader, opt_G, opt_D, device="cuda"): for low, high in loader: low, high = low.to(device), high.to(device) fake = G(low) pred_real = D(high) pred_fake = D(fake.detach()) real_label = torch.ones_like(pred_real, device=device) fake_label = torch.zeros_like(pred_fake, device=device) loss_D = (nn.functional.binary_cross_entropy_with_logits(pred_real, real_label) + nn.functional.binary_cross_entropy_with_logits(pred_fake, fake_label)) * 0.5 opt_D.zero_grad() loss_D.backward() opt_D.step() pred_fake = D(fake) loss_G = (nn.functional.binary_cross_entropy_with_logits(pred_fake, real_label) + 5.0 * nn.functional.l1_loss(fake, high) + 1.0 * decay_loss(fake, high)) opt_G.zero_grad() loss_G.backward() opt_G.step()参数说明：λ_l1=5.0是为了让生成矩阵和真实矩阵的绝对数值分布对齐，λ_decay=1.0是矫正距离衰减曲线。两个超参如果调大，生成结果会更保守，边界会变模糊；调小则容易生成一堆假信号。每次迭代先更新判别器，再更新生成器，这个顺序保证 D 不会太弱，对抗训练才有效。Notebook 里跑这段循环时，每 200 个 batch 打印一次 loss 数字即可，不要每个 batch 都 print，否则网页版的输出缓冲会越来越大，拖慢交互。
主机级别的超参可以按下面这张表起步：
| 超参数 | 取值 | 作用 | 调参方向 |
|---|---|---|---|
| λ_l1 | 5.0 | 保持整体强度 | 偏大则结果模糊，偏小则伪峰变多 |
| λ_decay | 1.0 | 约束距离衰减 | 偏大则对角线方向过度平滑 |
| G 学习率 | 2e-4 | 生成器更新速度 | 太高容易模式崩溃 |
| D 学习率 | 1e-4 | 判别器更新速度 | 太低则判别不出真假 |
4. Jupyter 环境下训练、可视化与排错
Notebook 和训练代码集成时，最大的敌人不是模型不收敛，而是环境不稳定。经验不多的时候，先用小尺寸patch_size=32、factor=2跑通全流程，确认种子、设备、可视化都正常，再换成正式参数量。
4.1 随机种子和硬件状态的统一管理
import random import numpy as np import torch def seed_everything(seed=42): random.seed(seed) np.random.seed(seed) torch.manual_seed(seed) torch.cuda.manual_seed_all(seed) seed_everything(42) device = "cuda" if torch.cuda.is_available() else "cpu" print("running on", device)把随机种子集中在函数里，课程设计里对比 baseline 和 GAN 增强效果时需要严格的可复现性。torch.cuda.manual_seed_all只对单机单卡有意义，多卡训练还要额外设置 DataLoader 的worker_init_fn，但一个 64x64 patch 的任务用不到多卡。
4.2 三个值得监控的训练指标
训练中间每 5 个 epoch 采样一组可视化，我会生成一张三合一图：低分辨率输入、GAN 输出、真实高分辨率。用 matplotlib 内联显示时注意三个图颜色映射和 vmin/vmax 保持一致，不然肉眼对比会失真。
import matplotlib.pyplot as plt def visualize_batch(low, fake, high, idx=0): fig, axes = plt.subplots(1, 3, figsize=(12, 4)) for ax, m, title in zip(axes, [low, fake, high], ["low input", "GAN output", "ground truth"]): im = ax.imshow(m[idx, 0].detach().cpu().numpy(), cmap="Reds", vmin=0, vmax=1) ax.set_title(title) ax.axis("off") plt.colorbar(im, ax=axes, fraction=0.04) plt.show()除了肉眼观察，距离衰减曲线是最适合课程设计的量化指标。它把接触强度按线性距离分组平均，直接反映 GAN 有没有保留长距离相互作用。
def contact_curve(mat): if torch.is_tensor(mat): mat = mat.squeeze().cpu().numpy() n = mat.shape[0] curve = np.zeros(n) for k in range(n): diag = np.diagonal(mat, offset=k) curve[k] = np.nanmean(diag) if len(diag) else 0.0 return curve把真实矩阵、低分辨率输入和 GAN 输出三条曲线画在同一张图里，看 GAN 的曲线是否向真实曲线靠拢。这里有一个很隐蔽的坑：生成的矩阵有时会在对角线方向出现断层，接触强度在某个位置突然降低，原因是 patch 边缘被卷积层抹平了，生成器只学会了局部连续性，没学会全局衰减。出现这个现象时，优先检查decay_loss的权重是否过小，以及切块的stride是否过小导致相邻 patch 彼此独立。
4.3 常见报错与处理方式
| 症状 | 可能原因 | 处理方式 |
|---|---|---|
| Jupyter Notebook 单元格执行代码没有任何反应 | 内核无响应，通常是显存不足或 numpy 版本冲突导致内核被杀 | 重启内核并看终端输出，先用小 batch 定位 |
| Jupyter Notebook 无法运行 | 端口被占用或旧版 notebook 与新版 Python 兼容问题 | 换 8890 等新端口，或pip install --upgrade notebook |
ImportError: DLL load failed while importing rpds | Windows 下 jsonschema 依赖的 rpds 二进制包损坏 | pip install --force-reinstall rpds-py jsonschema |
| CUDA out of memory | patch 太大或 batch_size 太大 | 把 patch_size 降回 64，或减小 batch_size |
| 生成结果像低分辨率输入的锐化版 | 判别器太弱，或 λ_l1 过大 | 先单独训练判别器 50 步，再把 λ_l1 降到 2.0 |
排错表的最后一行值得展开。很多同学看到生成结果“还行”就以为训练成功，其实生成器只是在学恒等映射。判断方法是把低分辨率输入和生成输出逐像素相减，如果差值非常小，说明生成器完全退化成插值器，对抗部分没有起作用。出现这种情况，把判别器的学习率从 1e-4 提到和生成器相同的 2e-4，通常能救回来。
如果你习惯在 nvim 里写代码，也可以起一个 Jupyter 内核后，通过jupyter-console或编辑器插件连到同一个内核上，交互式调试和绘图与网页版完全一致，不用在两个界面之间来回切换。
5. 期末大作业的验收清单与一个实用技巧
课程设计评阅最看重两点：任务定义是否合理，评价是否闭环。我建议报告里固定放一张对比表，把增强前后模型的量化指标放在一起，而不是只贴生成图。
| 评价维度 | 指标 | 数值参考 |
|---|---|---|
| 像素级相似度 | SSIM | 大于 0.90 算较好，低于 0.85 要检查归一化 |
| 线性相关性 | Pearson / Spearman | 大于 0.90，优先看 Spearman，对极值不敏感 |
| 距离衰减保持 | 曲线斜率和截距 | 对数坐标下拟合失败需要重建训练对 |
| 结构保持 | TAD 边界 F1 | 对真实矩阵和生成矩阵分别调用边界再对齐 |
SSIM 直接用skimage.metrics里的structural_similarity计算，先把两个矩阵归一化到 0-1，win_size设 7 或 11，避免窗口太小对对角线噪声敏感。
最后给一个验证增强效果的实用技巧：用绝缘分数（insulation score）检验生成矩阵有没有保留染色质边界结构。
def insulation(mat, window=10): n = mat.shape[0] score = np.zeros(n) for i in range(window, n - window): top_left = mat[i - window:i, i - window:i] bottom_right = mat[i + 1:i + 1 + window, i + 1:i + 1 + window] score[i] = top_left.mean() - bottom_right.mean() return score实现逻辑是沿对角线滑动一个窗口，对比左上方块和右下方块的接触强度差异。真实 Hi-C 矩阵在 TAD 边界两侧接触强度有明显跳变，绝缘分数的局部峰就对应边界位置。生成矩阵算出的峰如果和真实矩阵高度重合，说明增强没有破坏拓扑关联结构，这比单独看 SSIM 更有说服力，也正好覆盖老师对“方法对不对”的追问。注意这个绝缘分数是简化的一维版本，完整版一般会对两侧窗口分别归一化再取比值，但课程设计用这个简化版足够展示增强质量了。
跨染色体验证是最后一个技巧：用 chr7 训练的模型直接生成 chr19 的高分辨率矩阵，把绝缘分数和距离衰减曲线同样画出来，数值落差不超过 10% 才说明方法真正学到了普遍规律。模型文件命名建议是gan_hiC_chr7_f2_b32.pth，把染色体、factor、base 通道数嵌进文件名，训练日志里保留每个 epoch 的 loss 和指标记录。提交时把切块脚本、训练循环和读取代码放进同一个目录，评审直接双击 Notebook 就能从 .cool 文件出发重跑整个流程，这比任何文字描述都有说服力。
本文还有配套的精品资源，点击获取
