用于 RNAseq 比对的 Bash For Loop

Question

我的头脑根本不适用于 for 循环，因此我们将不胜感激。

背景：我正在尝试分析一些 RNAseq 数据，需要编写一个 for 循环来通过 STAR 读取我所有的双端 fastq 文件。

这是我现在的代码：

#! /bin/bash

#PBS -l walltime=24:00:00
#PBS -l nodes=1:ppn=12
#PBS -l mem=48

# ENVIRONMENT
module load gcc/6.2.0
module load STAR/2.6.1d
prj=/gpfs/data/elf-lab/RNA_Seq/RNA_trial/sab_sandbox/

# INPUTS
sample=SE-QS-19-FC01_S*

staridx=${prj}/reference/hg38_gencode28_STAR

fq1=${prj}/data/${sample}.R1_001.fastq
fq2=${prj}/data/${sample}.R2_001.fastq

rgline="ID:${sample}    PU:${sample}    SM:${sample}    PL:ILLUMINA LB:${sample}"

# OUTPUTS
outprefix=${prj}/alignment/${sample}.

# COMMAND
STAR \
  --runThreadN 12 \
  --genomeDir $staridx \
  --readFilesIn $fq1 $fq2 \
  --outSAMtype BAM Unsorted \
  --outSAMunmapped Within \
  --outFileNamePrefix $outprefix \
  --outSAMattrRGline $rgline \
  --outSAMattributes NH HI AS nM NM \
  --quantMode GeneCounts

这是我的文件的样子：

SE-QS-19-FC01_S33_R1_001.fastq.gz
SE-QS-19-FC01_S33_R2_001.fastq.gz
SE-QS-20-FC01_S34_R1_001.fastq.gz
SE-QS-20-FC01_S34_R2_001.fastq.gz

我想编写一个 for 循环，以便 fq1 和 fq2 将成为每次读取的每一对，但我不确定将 for 循环放置在何处以便可以在 STAR 命令中使用 fq1 和 fq2。提前谢谢你。

Answer 1

在使用数组之前，我已经遍历了 R1.fastq/R2.fastq 文件，假设您有相同数量的 R1/R2 文件，例如

#! /bin/bash

#PBS -l walltime=24:00:00
#PBS -l nodes=1:ppn=12
#PBS -l mem=48

# ENVIRONMENT
module load gcc/6.2.0
module load STAR/2.6.1d
prj=/gpfs/data/elf-lab/RNA_Seq/RNA_trial/sab_sandbox/

staridx=${prj}/reference/hg38_gencode28_STAR

## Make arrays named fq1 and fq2 ##
fq1=(${prj}/data/*.R1_001.fastq)
fq2=(${prj}/data/*.R2_001.fastq)

# COMMAND
for ((i=0;i<"${#fq1[@]}";i++)); do
  sample="${fq1[$i]%%_R*}"
  rgline="ID:${sample}    PU:${sample}    SM:${sample}    PL:ILLUMINA LB:${sample}"
  outprefix=${prj}/alignment/"${sample}"

  STAR \
  --runThreadN 12 \
  --genomeDir $staridx \
  --readFilesIn "${fq1[$i]}" "${fq2[$i]}" \
  --outSAMtype BAM Unsorted \
  --outSAMunmapped Within \
  --outFileNamePrefix $outprefix \
  --outSAMattrRGline $rgline \
  --outSAMattributes NH HI AS nM NM \
  --quantMode GeneCounts
done

如果您的文件已按照https://www.biostars.org/p/243683/ 进行 gzip 压缩（如您的示例中所示），您可能还需要添加 --readFilesCommand gunzip -c。

在这里进一步解释这个 c 风格的 for 循环：https://linuxize.com/post/bash-for-loop/#the-c-style-bash-for-loop

希望这会有所帮助 - 有多种方法可以解决此问题 - 如果您需要更多帮助，请在您的问题中添加更多详细信息，例如包括你用过但没用的代码。