如何从PLINK二进制BED/BIM/FAM文件中随机抽取3000个体子集

方法1：PLINK实现（推荐大文件场景，效率最高） #

PLINK是处理基因型二进制文件的原生工具，不需要将全量数据加载到内存，运行速度最快，适合十万级以上样本量的文件操作。

• 高版本PLINK（1.90b3.38及以上、PLINK2.0）可以直接用内置随机抽样参数一步完成：

plink --bfile 输入文件共同前缀 \
      --rand-sample 3000 \
      --seed 20240520 \
      --make-bed \
      --out 抽样后输出文件共同前缀

参数说明：

• --bfile：指定输入的BED/BIM/FAM三个文件的共同前缀
• --rand-sample 3000：无放回随机抽取3000个个体
• --seed：指定随机数种子，设置后可复现抽样结果，不需要可删掉
• --make-bed：指定输出格式为PLINK二进制BED/BIM/FAM格式
• 如果使用的PLINK版本较低不支持--rand-sample参数，可以先手动提取样本ID再做子集筛选：

# 从fam文件无放回随机抽3000个样本的FID和IID，存到keep文件
shuf -n 3000 输入文件共同前缀.fam | awk '{print $1,$2}' > keep_3000.txt
# 按keep列表提取样本
plink --bfile 输入文件共同前缀 \
      --keep keep_3000.txt \
      --make-bed \
      --out 抽样后输出文件共同前缀

方法2：R实现 #

适合已经在R环境中开展后续分析的场景，推荐用genio包操作PLINK二进制文件，接口简洁。

• 先安装依赖（未安装时执行）：

install.packages("genio")

• 抽样代码：

library(genio)
# 设定随机种子保证结果可复现
set.seed(20240520)
# 先读取fam文件获取全部样本信息，无需加载全量基因型，速度快
fam <- read_fam("输入文件共同前缀.fam")
# 随机抽取3000个样本的索引
keep_idx <- sample(x = nrow(fam), size = 3000, replace = FALSE)
# 读取BIM表和基因型矩阵
bim <- read_bim("输入文件共同前缀.bim")
geno <- read_bed("输入文件共同前缀.bed", n_ind = nrow(fam), n_snp = nrow(bim))
# 子集提取
fam_sub <- fam[keep_idx, ]
geno_sub <- geno[, keep_idx]
# 写出新的PLINK二进制文件，注意genio要求基因型矩阵为个体在行、SNP在列，因此做转置
write_plink(
  path = "抽样后输出文件共同前缀",
  X = t(geno_sub),
  bim = bim,
  fam = fam_sub
)

注意：如果总样本量超过10万，将全量基因型读入R会占用大量内存，这种情况优先选择PLINK方法。

方法3：Python实现 #

适合基于Python做下游分析的场景，使用pandas-plink库读写PLINK二进制文件。

• 先安装依赖：

pip install pandas pandas-plink numpy

• 抽样代码：

import numpy as np
import pandas as pd
from pandas_plink import read_plink, write_plink

# 设定随机种子保证结果可复现
np.random.seed(20240520)
# 读取PLINK文件，返回BIM表、FAM表、基因型矩阵（SNP在行，样本在列）
bim, fam, geno = read_plink("输入文件共同前缀")
# 随机抽取3000个样本的索引
keep_idx = np.random.choice(a = fam.shape[0], size = 3000, replace = False)
# 子集提取
fam_sub = fam.iloc[keep_idx, :]
geno_sub = geno[:, keep_idx]
# 写出新的PLINK二进制文件
write_plink(
  path = "抽样后输出文件共同前缀",
  bim = bim,
  fam = fam_sub,
  G = geno_sub
)

注意：超大规模文件读入Python时内存占用较高，优先选择PLINK方案提升效率。

结果校验 #

抽样完成后可直接检查输出的FAM文件行数确认样本量正确，Linux/macOS系统可直接运行命令校验：

wc -l 抽样后输出文件共同前缀.fam

返回结果为3000即代表抽样正确。

内容的提问来源于stack exchange，提问作者Yasaman Taher