1 叶绿体的组装

1.1 GetOrganelle的安装

此处默认Anaconda已经安装，那么使用下面的命令即可安装GetOrganelle：

conda create -n getorganelle python=3.7 #创建名叫getorganelle的环境
conda install -n getorganelle -c bioconda getorganelle #安装getorangelle

1.2 GetOrganelle的配置

安装好后需要配置参考基因组，GetOrganelle作者推荐以下命令：

get_organelle_config.py --add embplant_pt,embplant_mt 

然而该命令下载很慢，可以考虑到本地下载，手动导入，具体参考命令如下：

wget https://github.com/Kinggerm/GetOrganelleDB/releases/download/0.0.1/v0.0.1.tar.gz
tar -zxvf v0.0.1.tar.gz
get_organelle_config.py -a embplant_pt,embplant_mt --use-local ./0.0.1

该命令参考自：零基础教程 | 叶绿体基因组组装 - GetOrganelle

之后即可开始执行组装，但是可能会提示：

ERROR: Bowtie2 is not available!

实际上Bowtie2已经安装，但是似乎没有连接上，参考零基础教程 | 叶绿体基因组组装 - GetOrganelle的意见，直接更新bowtie2到最新版即可，实测也的确有效。

conda update bowtie2

1.3 使用GetOrganelle组装叶绿体

准备好自己的数据，一般为双末端的illumina数据，然后用以下命令运行：

get_organelle_from_reads.py -1 forward.fq -2 reverse.fq -o plastome_output 
    -R 15 -k 21,45,65,85,105 -F embplant_pt # 该代码来自getorganelle作者官网

-1 -2 输入双向测序数据，解压或者压缩文件均可
-o    指定输入文件的路径

-R    似乎是sample的循环迭代次数
-k    K-mer参数，不可超过reads的最大读长

-F    指定参考基因组，一般植物基因组为embplant_pt

开始运行后，可能会有一个SPAdes不会正常调用的报错：

......
2021-09-19 00:51:50,106 - ERROR: Error with running SPAdes: == Error ==  system call for: "['/home/zz/.conda/envs/getorganelle/share/spades-3.13.0-0/bin/spades-core', '/home/zz/genome/chloroplast/GetOrganelle-master/dinganensis/seed/embplant_pt.initial.fq.spades/K45/configs/config.info']" finished abnormally, err code: -11
......

该报错导致程序中断无法运行，经过google检索，发现似乎是个挺普遍的问题，参考Spades finished abnormally, err code: -11 · Issue #83 · Kinggerm/GetOrganelle · GitHub里面的建议，本地安装SPAdes的最新版本即可，因为 conda install spades -c bioconda 默认安装的不是最新版，奇怪的是bioconda源的确是有最新版的包的。

考虑到bioconda源中有最新版本的包，尝试直接在线安装指定版本：

conda install spades=3.15.3  

成功安装，getorganelle终于成功运行。

1.4 GetOrganelle组装叶绿体的参数优化

只要重测序的数据有一定的覆盖度（>2X），官方的命令一般即可组装成环，但是也常常出现不成环的情况，然后就需要对参数优化或者后续深入分析。

参数优化分为两个方向，一是增加reads的使用量和迭代运算量，二是使用自己的参考基因组：

# 增加数据使用量和迭代细节
get_organelle_from_reads.py -1 forward.fq -2 reverse.fq -o plastome_output 
    -R 30 -k 21,35,45,65,85,105,121,135,145  --reduce-reads-for-coverage 10000000 
    -F embplant_nr --overwrite
- R -k 提高迭代循环和k-mer粒度，也可设置更高的数值 
    --reduce-reads-for-coverage 10000000 个reads或者inf使用全部的reads
    这些设置可能使组装成环，但是运算量会成倍增加           

## 使用自己的参考基因组
get_organelle_from_reads.py -1 forward.fq -2 reverse.fq -o plastome_output 
    -R 15 -k 21,45,65,85,105 -F embplant_pt -s personal-reference.fasta
-s 指定自己的参考基因组，比如同一个物种、属的参考基因组，可以同时指定多个，放到同一个fasta文件即可

如果通过调整GetOrganelle参数仍然不能成环，可以使用Bandage来手动选择、删减，组合不同的config生成成环的组转结果，大致界面如下。具体可以参考此视频，或者Bilibili上的搬用视频，后续我们也可能考虑写个教程。

2. GetOrganelle的批量组装

如果批量测序了一批数据，逐个组装似乎有些麻烦，因此考虑批量组装，利用Linux的bash命令如下：

while read a; 
do get_organelle_from_reads.py -1 "$a"R1.fastq -2 "$a"R2.fastq -o "$a"-output 
        -R 15 -k 21,45,65,85,105 -F embplant_pt; 
done < list.txt
# 序列文件的前缀放到list.txt即可，会将list文件中的每一行作为a值，替换到do后面的命令中

3. GetOrganelle多个组装结果的合并

GetOrganelle的成环结果文件中通常同时包括两个构型的结果，但两种构型的顺序似乎是随机的。那么多个结果的合并就有些麻烦，需要比对后根据比对结果重选挑选另一个构型的结果。一般来说，不同构型结果的混合比对，速度相当缓慢。那么，使用以下的半成品小脚本可以用来合并多个GetOrganelle的结果，保证他们的构型一直，且同时可以检查不成环的结果。

# Author： Zhang Zhen
# Email： ficuszz@163.com
# Website: zhangzhen.plus
# Description： 这是一个半成品的脚本，需要自行调整一些诸如工作路径、文件输入路径、部分阈值之类的
#               参数。所以这里假设需要了解初步的R语言知识


library(ape)
library(strataG)

setwd("g:/Ubuntu_1804.2019.522.0_x64/rootfs/home/zz/byk/138-202209/") # 指定测序文件夹，该文件夹中应该只有多个样品测序的结果文件夹

file_path <- list.files()
sample_name <- sapply(file_path, function(x) substr(x, 1, nchar(x)-14)) # 此处的14为需要修剪掉的字符，具体随自己测序文件的命名规则而变
sample_number <- length(list.files()) 

# 将第一个样品的其中一条组装结果作为基准，将其余样品的同构型结果收集起来（请确保第一个样品组装成环）
j <- 1 # 如果需要选择另一个构型，请改为2
if(length(list.files(file_path[1], pattern = "fasta", full.names = T )) == 2 ) {
  temp0 <- read.FASTA(list.files(file_path[1], pattern = "fasta", full.names = T )[j])
  names(temp0) <- sample_name[1]
  write.FASTA(temp0, paste("c:/Users/qbycs/Desktop/byk/",sample_name[1], ".fasta", sep = ""))
} else {print("请确保第一个样品组装成环") }

#将其余样品的同构型组装结果筛选出来
checklist <- c()
for(i in 1:sample_number){
  if(length(list.files(file_path[i], pattern = "fasta", full.names = T )) == 2){ # 检查是否成环组装    
    temp <- read.FASTA(list.files(file_path[i], pattern = "fasta", full.names = T)[1])
    # 根据比对后遗传距离来判断构型
    if(dist.dna(mafft(c(temp0,temp))) < 0.01){ # 认为和基部遗传距离相差超过0.01为不同构型，然后自动选择另一条，可根据类群远近调整数值
      names(temp) <- sample_name[i] # 修改fasta的label
      write.FASTA(temp, paste("c:/Users/qbycs/Desktop/byk/",sample_name[i], ".fasta", sep = ""))
    } else {
      temp <- read.FASTA(list.files(file_path[i], pattern = "fasta", full.names = T )[2])
      names(temp) <- sample_name[i] # 修改fasta的label
      write.FASTA(temp, paste("c:/Users/qbycs/Desktop/byk/",sample_name[i], ".fasta", sep = ""))
    }
  } else {
    checklist <- c(checklist, sample_name[i])
    next;}    
}

# 输出问题组装样品名称
if(length(checklist) > 0){
  print(paste("请检查以下样品号的组装结果：", as.character(checklist)))
}