目录：
1. Data Input
2. Calculating Genotype Probabilities
3. Performing a genome scan
4. Performing a permutation test
5.
ok，进入正文。
1. Data Input
和 R/qtl 要求的数据格式不同，R/qtl2 要求每一个数据为一个单独的文件，但可以最后总结到一个 YAML 或 JSON 格式的文本里，一次性读入，还是挺方便的。
首先看一下这个 YAML 或 JSON 格式的文本文件，指定了群体类型，每一个文件的名字，以及基因型信息等。
其中，geno、pheno、gmap 三个文件是必须的，其它文件是可选的。
读入使用 read_cross2 函数。
read_cross2(file, quiet=TRUE)
# file 指定 YAML 文件的路径位置即可，quiet 指定是否输出读入过程。
以下分别是这 3 个文件的格式示例：
2. Calculating Genotype Probabilities
类似 qtl 包中的分析方法，QTL 分析首先需要计算已有标记之间的基因型概率。
使用的函数为 calc_genoprob() 。不过在使用之前，需要先在 markers 之间插入 pseudomarkers，以确定要计算基因型可能性的位置。在遗传图谱中插入 pseudomarkers 的函数为 insert_pseudomarkers() 。
举例示范一下：
library(qtl2)
## 读入包里自带的示例文件
iron <- read_cross2(file = system.file("extdata", "iron.zip", package="qtl2") )
## 在遗传图谱内插入 pseudomarkers
map <- insert_pseudomarkers(map=iron$gmap, step=1)
# step 指的是 pseudomarkers 和标记之间的举例为 1 cM
## 计算 pseudomarkers 的基因型概率
pr <- calc_genoprob(cross=iron, map=map, error_prob=0.002, cores=4)
# error_prob 指定假设的基因分型错误概率，默认为0.0001；cores 是在集群上指定运算使用的内核数，parallel::detectCores() 查看可用的内核数
这个过程就到此为止了，接下来通过几张图看一下每一步的效果。
summary(iron)
head(iron$gmap)
map <- insert_pseudomarkers(map=iron$gmap, step=1)
head(map, n=1) # 查看1号染色体上的标记
pr <- calc_genoprob(cross=iron, map=map, error_prob=0.002)
(pr$`19`)[1:3,,"c19.loc4"] # 查看19号染色体上该标记处前3个个体的基因型概率
此外，如果需要使用加性模型进行基因组扫描，需要将基因型概率转换为等位基因概率。
apr <- genoprob_to_alleleprob(probs=pr)
3. Performing a genome scan
原理就不在这块详细说了，有兴趣的可以去看 R/qtl 定位分析（三）Single-QTL analysis 。
以 Haley-Knott regression 方法为例：
Xcovar <- get_x_covar(iron)
out <- scan1(genoprobs = pr, pheno = iron$pheno, Xcovar=Xcovar, cores=4)
输出是 LOD 得分矩阵，位置×表型。
head(out, n=10)
plot_scan1(out, map = map, lodcolumn = "liver")
4. Performing a permutation test
虽然可以看到部分染色体上存在峰，但是否显著？对结果如何评估？
类似于 qtl 包中的操作，需要进行排列检验，打乱基因型和表型，去除二者之间的关联，然后进行计算全基因组的最大 LOD 评分作为阈值，从而确定结果是否是随机出现的。需要用到的函数为 scan1perm() 。
operm <- scan1perm(genoprobs = pr, pheno = iron$pheno, Xcovar = Xcovar, n_perm = 1000)
# n_perm 指定排列检验的重复次数
hist(operm[,'liver'], breaks = 50, xlab = "LOD", main = "LOD scores for liver scan with threshold in red")abline(v = summary(operm)[,'liver'], col = 'red', lwd = 2)
summary(operm, alpha=c(0.1, 0.05))
如果要对 X 染色体单独计算阈值，可使用下面的命令：
operm2 <- scan1perm(pr, iron$pheno, Xcovar=Xcovar, n_perm=1000, perm_Xsp=TRUE, chr_lengths=chr_lengths(map))
summary(operm2, alpha=c(0.2, 0.05))
5. Finding LOD peaks
这一步就要根据排列检验得到的阈值，以及计算结果，确定目标的 QTL 区间了。首先要确定 LOD peaks 的位置，再计算其可靠的区间范围。
operm <- scan1perm(pr, iron$pheno, Xcovar=Xcovar, n_perm=1000)
thr <- summary(operm)
find_peaks(scan1_output = out, map = map, threshold = thr, prob = 0.95, expand2markers = FALSE)
# prob 确定贝叶斯置信区间，0.95 指的是有 95% 的可能性落在该区间
# expand2markers 指的是是否扩展至相邻标记
该函数也可以确定一条染色体上的多个峰，不过需要 peakdrop 指定两个相邻峰值之间的最低值必须下降的量。
find_peaks(scan1_output = out, map = map, threshold = thr, peakdrop = 1.8, prob = 0.95, expand2markers = FALSE)
