step3 megahit组装

困难重重的拼接，按照软件教程，一直报错

报错原因

随后咨询了朋友，seqtk排除文件损坏，排除内存不足，最后发现是环境不匹配。用下面代码就成功了

singularity exec docker://quay.io/biocontainers/megahit:1.2.9--h5b5514e_3 \
    megahit -1 SRR22201374_NDF20_paired_1.fastq \
            -2 SRR22201374_NDF20_paired_2.fastq \
            -o megahit_container_output \
            -t 4

随后使用脚本

因为上一步每一个样本在一个文件夹，我只需要去掉宿主基因组后的fastq文件
于是建立了assembly_input文件夹把文件软链接过来
脚本如下

#!/bin/bash

#############
##1、软链接数据
kneaddata="/home/shpc_101389/wang_proj/yang_data/kneaddata"
assembly_input="/home/shpc_101389/wang_proj/yang_data/assembly_input"
mkdir -p "$assembly_input"

for d in "$kneaddata"/*/; do
  sample=$(basename "$d")

 r1=( "$d"/*paired_1.fastq )
  r2=( "$d"/*paired_2.fastq )

  if (( ${#r1[@]} == 0 || ${#r2[@]} == 0 )); then
    echo "跳过 $sample：找不到 *paired_1.fastq / *paired_2.fastq"
    continue
  fi

 ln -sf "${r1[0]}" "$assembly_input/${sample}_R1.fastq"
  ln -sf "${r2[0]}" "$assembly_input/${sample}_R2.fastq"
  echo "已链接：$sample"
done



echo "完成。软链接目录：$assembly_input"


OUT_ROOT="/home/shpc_101389/wang_proj/yang_data/assembly_output/megahit_container_output"
IMAGE="docker://quay.io/biocontainers/megahit:1.2.9--h5b5514e_3"
THREADS=10
MIN_CONTIG=500

mkdir -p "$OUT_ROOT"
cd "$assembly_input"

for r1 in *R1.fastq *R1.fastq.gz; do
  r2="${r1/R1.fastq/R2.fastq}"
  r2="${r2/R1.fastq.gz/R2.fastq.gz}"
  [[ -f "$r2" ]] || { echo "skip (no R2): $r1"; continue; }

#样本名与输出目录
  sample="${r1%R1.fastq}"; sample="${sample%R1.fastq.gz}"
  out="$OUT_ROOT/${sample%[_\.]}"

 echo "[$(date '+%F %T')] MEGAHIT start: $sample"

 singularity exec -B "$assembly_input:$assembly_input" -B "$OUT_ROOT:$OUT_ROOT" \
    "$IMAGE" megahit \
      -t "$THREADS" \
      -1 "$assembly_input/$r1" \
      -2 "$assembly_input/$r2" \
      -o "$out" \
      --min-contig-len "$MIN_CONTIG" \
    2>&1 | tee "$OUT_ROOT/${sample}.log"

  echo "[$(date '+%F %T')] MEGAHIT done: $sample"
  echo "结果目录：$out"
  echo "contigs：$out/final.contigs.fa"
done

echo "MEGAHIT组装完成！结果保存在：$OUT_ROOT"

最终每个样本经过拼接后生成文件如图

拼接结果.png

安装quast评估质量

一直安不上重新建个环境

conda create -n quast -c conda-forge -c bioconda --strict-channel-priority quast
conda activate quast#以后使用记得激活这步环境

提取所有fa文件
quast.py *.fa