Multiplexed CRISPR-based target-enriched next-generation sequencing for detecting antibiotic resistance genes in environmental samples V.2

import re

file=open(r"INPUT_FASTA_FILE_PATH_HERE")
raw_content=file.readlines()
file.close()

output_content=[]

for i in range(0,len(raw_content)):
    if i%2==0:
        aro_number=re.findall(r'ARO:[0-9]+',raw_content[i])[0]
        output_content.append('>'+aro_number.split(':')[1])
    if i%2!=0:
        output_content.append(raw_content[i].strip('\n'))

output_file=open(r" OUTPUT_TRIMMED_FASTA_FILE_PATH_HERE ",'w')

for i in range(0,len(output_content)):
    output_file.write(output_content[i]+'\n')

output_file.close()

split -b 10m INPUT_FILE_NAME OUTPUT_FILE_PATH_AND_PREFIX

ls | while read i; do mv ${i} ${i}.txt

import argparse
import os

parser = argparse.ArgumentParser()
parser.add_argument('-i',dest='input',type=str,required=True,help='Define input txt folder')
parser.add_argument('-o',dest='output',type=str,required=True,help="Define output path for a combined txt")
args=parser.parse_args()

full_gRNA_list=[]
dir_path=args.input.strip("'")
dir_list=os.listdir(dir_path)
for file_name in dir_list:
    print('Processing '+file_name+' ......')
    file=open(dir_path+'/'+file_name)
    raw_content=file.readlines()
    file.close()
    for i in range(0,len(raw_content)):
        raw_content[i]=raw_content[i].strip('\t\n\r')
    for i in range(0,len(raw_content),3):
        full_gRNA_list.append(raw_content[i])

full_gRNA_list=list(set(full_gRNA_list))

output_file=open(args.output.strip("'"),"w")
for i in range(0,len(full_gRNA_list)):
    output_file.write(full_gRNA_list[i]+'\n')
output_file.close()

import argparse

parser = argparse.ArgumentParser()
parser.add_argument('-i',dest='input',type=str,required=True,help='Define input txt folder')
parser.add_argument('-o',dest='output',type=str,required=True,help="Define output path for a combined txt")
args=parser.parse_args()

file=open(args.input.strip("'"))
gRNA_list=file.readlines()
file.close()

output_list=sorted(list(set(gRNA_list)))
output_file=open(args.output.strip("'"),'w')
for i in range(0,len(output_list)):
    output_file.write(output_list[i])
output_file.close()

java -jar PATH_TO/trimmomatic-0.39.jar SE --threads n input_read1_or_2.fastq.gz output_UMI.fastq.gz HEADCROP:8

cat file.fastq | paste - - - - | sort -k1,1 -t " " | tr "\t" "\n" > file_sorted.fastq

umi_tools extract --extract-method=string --bc-pattern=NNNNNNNNN -I [UMI_index_fastq_file] -S [UMI_index_fastq_file_out_nothing_left] --read2-in=[read_fastq_file_to_add_UMI] --read2-out=[output_read_fastq_file_with_UMI_added]

PATH_TO_PRICE_FOLDER/Price/PriceSeqFilter -fp R1.fastq R2.fastq -rqf 85 0.98 -rnf 90 -op R1_filtered.fastq R2_filtered.fastq

bowtie2-build reference.fasta reference

bowtie2 -x reference -1 R1_with_UMI.fastq -2 R2_with_UMI.fastq -a --very-sensitive --threads n | samtools view -bSF4 - | samtools sort - -o sorted_output.bam -@ n

samtools index sorted_output.bam -@ n   

umi_tools dedup -I input.bam --output-stats=deduplicated -S deduplicated.bam

samtools sort deduplicated.bam -o deduplicated_sorted.bam -@ n
samtools index deduplicated_sorted.bam -@ n  

samtools coverage --ff 0 -o coverage_out.tsv deduplicated_sorted.bam

import pandas as pd
import os
import argparse

parser = argparse.ArgumentParser()
parser.add_argument('-i',dest='input',type=str,required=True,help='Define input tsv folder')
parser.add_argument('-r',dest='reference',type=str,required=True,help='Define input reference file')
parser.add_argument('-o',dest='output',type=str,required=True,help="Define output path for a combined excel with ARG names and drug classes")
args=parser.parse_args()

ARO_table=pd.read_csv(args.reference,sep='\t')

dir_path=args.input.strip("'")
dir_list=os.listdir(dir_path)
dir_list=sorted(dir_list)
for file_name in dir_list:
    if file_name.endswith('.tsv'):
        print(file_name)
        prefix=file_name.split('.tsv')[0]        
        mapping_table=pd.read_csv(dir_path+'/'+file_name,sep='\t')
        mapping_table['ARG_name']=''
        mapping_table['ARG_family']=''
        mapping_table['Drug_class']=''
        
        for i in range(0,len(mapping_table)):
            if mapping_table.loc[i,'numreads']==0:
                continue
            else:
                for j in range(0,len(ARO_table)):
                    if str(mapping_table.loc[i,'#rname']).strip('\t ') in ARO_table.loc[j,'ARO Accession']:
                        mapping_table.loc[i,'ARG_name']=ARO_table.loc[j,'CARD Short Name']
                        mapping_table.loc[i,'Drug_class']=ARO_table.loc[j,'Drug Class']
                        mapping_table.loc[i,'ARG_family']=ARO_table.loc[j,'AMR Gene Family']
                        break
        
        output_table=pd.DataFrame(columns=mapping_table.columns)
        for i in range(0,len(mapping_table)):
            if mapping_table.loc[i,'coverage']>=75:
                output_table=output_table.append(mapping_table.loc[i],ignore_index=True)
             
    if os.path.exists(args.output):
        with pd.ExcelWriter(args.output, mode='a') as writer:  
            output_table.to_excel(writer, sheet_name=prefix) 
    else:
        output_table.to_excel(args.output,sheet_name=prefix)