Skip to content
Varsha V edited this page Jan 16, 2025 · 5 revisions

Welcome to the spore-to-VO wiki page. Here, you will find code snippets used for characterizing 14 novel species isolated from the NASA cleanroom Spacecraft Assembly Facility (SAF) at the Jet Propulsion Laboratory (JPL), California, during the Mars 2020 mission spacecraft assembly.

1. Genome Assembly

2. Overall Genome Related Indices (OGRI) Calculation

We retrieved all validly published representative genomes of the identified genera from NCBI using the command-line tool BIT (https://github.com/AstrobioMike/bit) along with NCBI Entrez-Direct.

# installation
conda install -c bioconda entrez-direct 
conda install -c conda-forge -c bioconda -c defaults -c astrobiomike bit 

# esearch to get the list of representative species accession IDs
esearch -query '<Agrococcus>[ORGN] AND "representative genome"[filter] AND all[filter] NOT anomalous[filter]' -db assembly | esummary | xtract - 
pattern DocumentSummary -def "NA" -element AssemblyAccession > <Agrococcus>_repr.txt

# bit to fetch the genomes
bit-dl-ncbi-assemblies -w <Agrococcus>_repr.txt -f fasta -j 100

# unzipping
gunzip *.gz

We determined the Average Nucleotide Identity (ANI) of the novel species in comparison to other type strains using FastANI version 1.34 (https://github.com/ParBLiSS/FastANI).

# get all the paths of the genomes
find ${genome_dir} -iname *.fasta > Novel_Genomes_Paths.txt

# perform all vs all ANI with the genomes
fastANI --rl Novel_Genomes_Paths.txt --ql Novel_Genomes_Paths.txt --threads 50 --matrix -o Novel_Genomes_ANI.tsv

Additionally, Average Amino Acid Identity (AAI) values were calculated using the aai.rb function from the Enveomics Collection toolbox (https://github.com/lmrodriguezr/enveomics).

aai.rb --seq1 <genome> --seq2 <ref_genome> --threads 50 --res AAI_output.txt

We annotated 14 novel species utilizing Prokka, a widely used command-line annotation tool (https://github.com/tseemann/prokka).

# prokka command
# From the directory of the .fasta files
for file in *.fasta; do prokka --outdir ../Prokka_Annotation/"$file" --prefix "$file" --cpus 40 "$file"; done

From the annotated genomes, we extracted the 16S rRNA and DNA gyrase subunit-B (gyrB) genes for the novel species and its closest type strain. Sequence identity was then assessed using Blast version 2.13.0.

# Extracting DNA gyrase subunit B
pattern_gyrb="DNA gyrase subunit B"

echo "Getting gyrB sequences..."

# Set Internal Field Separator (IFS) to handle newline properly in filenames
IFS=$'\n'

# Loop through all .ffn files in the input folder
for i in $(find "$INPUT_FOLDER" -iname "*.ffn"); do
  # Sanitize filename by replacing spaces with underscores and removing the extension
  filename=$(basename "$i" | sed -e 's/ /_/g' | cut -d . -f 1-2)
  echo "Getting gyrB sequence for: $i"

  # Use awk to extract the longest sequence that matches the pattern
  awk -v pattern="$pattern_gyrb" '
   BEGIN { RS=">"; FS="\n" }
   tolower($1) ~ tolower(pattern) {
    sequence = $0
    gsub(/\n/, "", sequence)  # Remove newline characters
    length_seq = length(sequence)
  
    # Update max_length and max_sequence if the current sequence is longer
    if (length_seq > max_length) {
      max_length = length_seq
     max_sequence = ">" $0
    }
   }
   END { if (max_sequence) printf "%s", max_sequence }' "$i" > "$OUTPUT_FOLDER/${filename}_gyrB.fasta"
done

echo "gyrB sequence extraction complete."

3. WGS-based phylogeny

We used GToTree version 1.8.2 (https://github.com/AstrobioMike/GToTree) to prepare phylogenetic trees based on single-copy core genes and subsequently performed 1,000 ultrafast bootstrap replicates with IQ-TREE2 version 2.2.0.3 (https://github.com/iqtree/iqtree2).

# GToTree
GToTree -a ${genus}_non_novel_accession.txt -f ${genus}_paths.txt -H $scg_set -t -L Species,Strain -T IQ-TREE -j 50 -m 
GToTree_Mapping_Files/${genus}_GToTree_Mapping.tsv -o $genus

# IQTREE-2
iqtree2 -s $(genus}/Aligned_SCGs_mod_names.faa -spp ${genus}/run_files/Partitions.txt -m txt -m MFP -bb 1000 -nt 40 -pre ${genus}_IQ_Tree

4. Genome Characterisation

The Cluster of Orthologous Genes (COGs) were identified using the genomes annotated with Prokka. This was achieved with the Python package cogclassifier (https://pypi.org/project/cogclassifier/).

# cogclassfier command
# Once Prokka is over, utilise .faa files for COG identification
for file in ./*/*.faa; do cp "$file" ./../COGs/; done

# From the directory of the .faa files
for files in *.faa; do COGclassifier -i "$files" -o ./COG_Output/"$files"; done

Subsequently, we predicted, annotated, and analyzed the secondary metabolite biosynthesis gene clusters in these novel strains using the command-line tool antiSMASH version 7.0 (https://antismash.secondarymetabolites.org/#!/about).

# antismash command
# From the directory of the .fasta files
for file in *.fasta; do antismash "$file" --taxon bacteria --cpus 20 --genefinding-tool prodigal --cb-knownclusters --output-dir ./../antiSMASH/"$file" --cc-mibig --fullhmmer; done

Additionally, we used the command-line tool Resistance Gene Identifier (RGI) (https://github.com/arpcard/rgi) to predict the antibiotic resistomes in these novel genomes. The identified antibiotic resistance genes (ARGs) were compared against the Comprehensive Antibiotic Resistance Database (CARD) (https://card.mcmaster.ca/).

# rgi command
# From the directory of the .fasta files
for file in *.fasta; do rgi main -i "$file" -o ./RGI_CARD/"$file" -t contig -a BLAST -n 50; done

5. Mapping metagenomics samples to novel isolates

We utilized fastp version 0.22.0 (https://github.com/OpenGene/fastp) to filter the raw shotgun sequences collected from NASA cleanrooms.

fastp.0.23.4 --in1 ${metagenomeFolder}/$read1 --out1 ${read1}_output.fastq --in2 ${metagenomeFolder}/{read2} --out2 ${read2}_output.fastq -w 16 --json ${metagenomeFolder}_fastp.json --html ${metagenomeFolder}_fastp.html

Next, we used MetaCompass version 2.0 (https://github.com/marbl/MetaCompass) to align the filtered reads to the newly identified genomes and estimate their abundance in NASA cleanrooms based on the mapped reads.

python3 ../../MetaCompass/go_metacompass.py -r $genome -1 "$read1" -2 "$read2" -l 150 -t 70 -y 50 -o Metagenomic_Mapping/${name}_${i}

Clone this wiki locally