Use miRDeep2 (Friedländer et al. 2011) to identify potential miRNAs using P.evermanni sRNAseq reads.
Inputs:
Requires collapsed reads (i.e. concatenated, unique reads) in FastA format. See 10.2-Peve-sRNAseq-BLASTn-31bp-fastp-merged.Rmd for code.
Genome FastA. See 07-Peve-sRNAseq-MirMachine.Rmd for download info if needed.
miRBase mature miRNAs FastA.
Utilizes a modified version, which includes cnidarian miRNA culled from literature by Jill Ahsley.
Outputs:
Primary outputs are a result table in BED, CSV (tab-delimited), and HTML formats.
Due to the nature of mirDeep2’s naming, trying to use variable names is challenging. As such, the chunks processing those files will require manual intervention to identify and provide the output filename(s); they are not handled at in the .bashvars file at the top of this script.
Other output files are too large for GitHub (and some (all?) are not needed for this analysis). Please find a full backup here:
This allows usage of Bash variables across R Markdown chunks.
{
echo "#### Assign Variables ####"
echo ""
echo "# Trimmed FastQ naming pattern"
echo "export trimmed_fastqs_pattern='*-31bp-merged.fq.gz'"
echo ""
echo "# Data directories"
echo 'export deep_dive_dir=/home/shared/8TB_HDD_02/shedurkin/deep-dive'
echo 'export deep_dive_data_dir="${deep_dive_dir}/data"'
echo 'export output_dir_top=${deep_dive_dir}/E-Peve/output/11.1-Peve-sRNAseq-miRdeep2-31bp-fastp-merged-cnidarian_miRBase'
echo 'export genome_fasta_dir=${deep_dive_dir}/E-Peve/data'
echo 'export trimmed_fastqs_dir="${deep_dive_dir}/E-Peve/output/06.2-Peve-sRNAseq-trimming-31bp-fastp-merged/trimmed-reads"'
echo 'export collapsed_reads_dir="${deep_dive_dir}/E-Peve/output/10.2-Peve-sRNAseq-BLASTn-31bp-fastp-merged"'
echo ""
echo "# Input/Output files"
echo 'export collapsed_reads_fasta="collapsed-reads-all.fasta"'
echo 'export collapsed_reads_mirdeep2="collapsed-reads-all-mirdeep2.fasta"'
echo 'export collapsed_reads_mirdeep2_over17bp="collapsed-reads-over17bp-mirdeep2.fasta"'
echo 'export concatenated_trimmed_reads_fastq="concatenated-trimmed-reads-all.fastq.gz"'
echo 'export genome_fasta_name="Porites_evermanni_v1.fa"'
echo 'export genome_fasta_no_spaces="Porites_evermanni_v1_genomic-no_spaces.fna"'
echo 'export mirdeep2_mapping_file="Peve-mirdeep2-mapping.arf"'
echo 'export mirbase_mature_fasta_name="cnidarian-mirbase-mature-v22.1.fasta"'
echo 'export mirbase_mature_fasta_no_spaces="cnidarian-mirbase-mature-v22.1-no_spaces.fa"'
echo ""
echo "# Paths to programs"
echo 'export mirdeep2_mapper="mapper.pl"'
echo 'export mirdeep2="miRDeep2.pl"'
echo 'export bowtie_build="/home/shared/bowtie-1.3.1-linux-x86_64/bowtie-build"'
echo ""
echo "# Set number of CPUs to use"
echo 'export threads=46'
echo ""
echo "# Initialize arrays"
echo 'export trimmed_fastqs_array=()'
} > .bashvars
cat .bashvars#### Assign Variables ####
# Trimmed FastQ naming pattern
export trimmed_fastqs_pattern='*-31bp-merged.fq.gz'
# Data directories
export deep_dive_dir=/home/shared/8TB_HDD_02/shedurkin/deep-dive
export deep_dive_data_dir="${deep_dive_dir}/data"
export output_dir_top=${deep_dive_dir}/E-Peve/output/11.1-Peve-sRNAseq-miRdeep2-31bp-fastp-merged-cnidarian_miRBase
export genome_fasta_dir=${deep_dive_dir}/E-Peve/data
export trimmed_fastqs_dir="${deep_dive_dir}/E-Peve/output/06.2-Peve-sRNAseq-trimming-31bp-fastp-merged/trimmed-reads"
export collapsed_reads_dir="${deep_dive_dir}/E-Peve/output/10.2-Peve-sRNAseq-BLASTn-31bp-fastp-merged"
# Input/Output files
export collapsed_reads_fasta="collapsed-reads-all.fasta"
export collapsed_reads_mirdeep2="collapsed-reads-all-mirdeep2.fasta"
export collapsed_reads_mirdeep2_over17bp="collapsed-reads-over17bp-mirdeep2.fasta"
export concatenated_trimmed_reads_fastq="concatenated-trimmed-reads-all.fastq.gz"
export genome_fasta_name="Porites_evermanni_v1.fa"
export genome_fasta_no_spaces="Porites_evermanni_v1_genomic-no_spaces.fna"
export mirdeep2_mapping_file="Peve-mirdeep2-mapping.arf"
export mirbase_mature_fasta_name="cnidarian-mirbase-mature-v22.1.fasta"
export mirbase_mature_fasta_no_spaces="cnidarian-mirbase-mature-v22.1-no_spaces.fa"
# Paths to programs
export mirdeep2_mapper="mapper.pl"
export mirdeep2="miRDeep2.pl"
export bowtie_build="/home/shared/bowtie-1.3.1-linux-x86_64/bowtie-build"
# Set number of CPUs to use
export threads=46
# Initialize arrays
export trimmed_fastqs_array=()Per miRDeep2 documentation:
The readID must end with _xNumber and is not allowed to contain whitespaces. has to have the format name_uniqueNumber_xnumber
# Load bash variables into memory
source .bashvars
sed '/^>/ s/-/_x/g' "${collapsed_reads_dir}/${collapsed_reads_fasta}" \
| sed '/^>/ s/>/>seq_/' \
> "${output_dir_top}/${collapsed_reads_mirdeep2}"
grep "^>" ${collapsed_reads_dir}/${collapsed_reads_fasta} \
| head
echo ""
echo "--------------------------------------------------"
echo ""
grep "^>" ${output_dir_top}/${collapsed_reads_mirdeep2} \
| head>1-2039267
>2-1218024
>3-551889
>4-472935
>5-375248
>6-292862
>7-292012
>8-258161
>9-226514
>10-214759
--------------------------------------------------
>seq_1_x2039267
>seq_2_x1218024
>seq_3_x551889
>seq_4_x472935
>seq_5_x375248
>seq_6_x292862
>seq_7_x292012
>seq_8_x258161
>seq_9_x226514
>seq_10_x214759mirDeep2 only accepts sequences at least 17bp long
# Load bash variables into memory
source .bashvars
# Only do this if you haven't already filtered by size
if [ ! -e "${output_dir_top}/${collapsed_reads_mirdeep2_over17bp}" ]; then
# Loop through the file
while IFS= read -r line || [[ -n "$line" ]]; do
# echo "checking line $line"
# Check if the line contains a sequence ID
if [[ "$line" == ">"* ]]; then
sequence_id="$line"
read -r sequence_line
# Check if the sequence length is greater than or equal to 17
if [ ${#sequence_line} -ge 17 ]; then
#echo "Writing to filtered_sequences_long.txt: $sequence_id"
echo "$sequence_id" >> "${output_dir_top}/${collapsed_reads_mirdeep2_over17bp}"
echo "$sequence_line" >> "${output_dir_top}/${collapsed_reads_mirdeep2_over17bp}"
else
#echo "Writing to filtered_sequences_short.txt: $sequence_id"
echo "$sequence_id" >> "${output_dir_top}/less_than_17bp.fasta"
echo "$sequence_line" >> "${output_dir_top}/less_than_17bp.fasta"
fi
fi
done < "${output_dir_top}/${collapsed_reads_mirdeep2}"
fi
wc -l ${output_dir_top}/${collapsed_reads_mirdeep2}
wc -l ${output_dir_top}/${collapsed_reads_mirdeep2_over17bp}
head ${output_dir_top}/less_than_17bp.fastamiRDeep2 can’t process genome FastAs with spaces in the description lines.
So, I’m replacing spaces with underscores.
And, for aesthetics, I’m also removing commas.
# Load bash variables into memory
source .bashvars
sed '/^>/ s/ /_/g' "${genome_fasta_dir}/${genome_fasta_name}" \
| sed '/^>/ s/,//g' \
> "${genome_fasta_dir}/${genome_fasta_no_spaces}"
grep "^>" ${genome_fasta_dir}/${genome_fasta_name} \
| head
echo ""
echo "--------------------------------------------------"
echo ""
grep "^>" ${genome_fasta_dir}/${genome_fasta_no_spaces} \
| head>Porites_evermani_scaffold_1
>Porites_evermani_scaffold_2
>Porites_evermani_scaffold_3
>Porites_evermani_scaffold_4
>Porites_evermani_scaffold_5
>Porites_evermani_scaffold_6
>Porites_evermani_scaffold_7
>Porites_evermani_scaffold_8
>Porites_evermani_scaffold_9
>Porites_evermani_scaffold_10
--------------------------------------------------
>Porites_evermani_scaffold_1
>Porites_evermani_scaffold_2
>Porites_evermani_scaffold_3
>Porites_evermani_scaffold_4
>Porites_evermani_scaffold_5
>Porites_evermani_scaffold_6
>Porites_evermani_scaffold_7
>Porites_evermani_scaffold_8
>Porites_evermani_scaffold_9
>Porites_evermani_scaffold_10# Load bash variables into memory
source .bashvars
sed '/^>/ s/ /_/g' "${deep_dive_data_dir}/${mirbase_mature_fasta_name}" \
| sed '/^>/ s/,//g' \
> "${deep_dive_data_dir}/${mirbase_mature_fasta_no_spaces}"
grep "^>" ${deep_dive_data_dir}/${mirbase_mature_fasta_name} \
| head
echo ""
echo "--------------------------------------------------"
echo ""
grep "^>" ${deep_dive_data_dir}/${mirbase_mature_fasta_no_spaces} \
| head>cel-let-7-5p MIMAT0000001 Caenorhabditis elegans let-7-5p
>cel-let-7-3p MIMAT0015091 Caenorhabditis elegans let-7-3p
>cel-lin-4-5p MIMAT0000002 Caenorhabditis elegans lin-4-5p
>cel-lin-4-3p MIMAT0015092 Caenorhabditis elegans lin-4-3p
>cel-miR-1-5p MIMAT0020301 Caenorhabditis elegans miR-1-5p
>cel-miR-1-3p MIMAT0000003 Caenorhabditis elegans miR-1-3p
>cel-miR-2-5p MIMAT0020302 Caenorhabditis elegans miR-2-5p
>cel-miR-2-3p MIMAT0000004 Caenorhabditis elegans miR-2-3p
>cel-miR-34-5p MIMAT0000005 Caenorhabditis elegans miR-34-5p
>cel-miR-34-3p MIMAT0015093 Caenorhabditis elegans miR-34-3p
--------------------------------------------------
>cel-let-7-5p_MIMAT0000001_Caenorhabditis_elegans_let-7-5p
>cel-let-7-3p_MIMAT0015091_Caenorhabditis_elegans_let-7-3p
>cel-lin-4-5p_MIMAT0000002_Caenorhabditis_elegans_lin-4-5p
>cel-lin-4-3p_MIMAT0015092_Caenorhabditis_elegans_lin-4-3p
>cel-miR-1-5p_MIMAT0020301_Caenorhabditis_elegans_miR-1-5p
>cel-miR-1-3p_MIMAT0000003_Caenorhabditis_elegans_miR-1-3p
>cel-miR-2-5p_MIMAT0020302_Caenorhabditis_elegans_miR-2-5p
>cel-miR-2-3p_MIMAT0000004_Caenorhabditis_elegans_miR-2-3p
>cel-miR-34-5p_MIMAT0000005_Caenorhabditis_elegans_miR-34-5p
>cel-miR-34-3p_MIMAT0015093_Caenorhabditis_elegans_miR-34-3pmiRDeep2 requires a Bowtie v1 genome index - cannot use Bowtie2 index
# Load bash variables into memory
source .bashvars
# Check for existence of genome index first
if [ ! -f "${genome_fasta_no_spaces%.*}.ebwt" ]; then
${bowtie_build} \
${genome_fasta_dir}/${genome_fasta_no_spaces} \
${genome_fasta_dir}/${genome_fasta_no_spaces%.*} \
--threads ${threads} \
--quiet
fiRequires genome to be previously indexed with Bowtie.
Additionally, requires user to enter path to their mirdeep2 directory as well
as their perl5 installation.
# Load bash variables into memory
source .bashvars
# Append miRDeep2 to system PATH and set PERL5LIB
export PATH=$PATH:/home/shared/mirdeep2/bin
export PERL5LIB=$PERL5LIB:/home/shared/mirdeep2/lib/perl5
# Run miRDeep2 mapping
time \
${mirdeep2_mapper} \
${output_dir_top}/${collapsed_reads_mirdeep2_over17bp} \
-c \
-p ${genome_fasta_dir}/${genome_fasta_no_spaces%.*} \
-t ${output_dir_top}/${mirdeep2_mapping_file} \
-o ${threads}Recommendation is to use the closest related species in the miRDeep2 options, even if the species isn’t very closely related. The documentation indicates that miRDeep2 is always more accurate when at least a species is provided.
The options provided to the command are as follows:
none: Known miRNAs of the species being analyzed.-t S.pupuratus: Related species.none: Known miRNA precursors in this species.-P: Specifies miRBase version > 18.-v: Remove temporary files after completion.-g -1: Number of precursors to anlayze. A setting of -1 will analyze all. Default is 50,000
I set this to -1 after multiple attempts to run using the default kept failing.NOTE: This will take an extremely long time to run (days). Could possible by shortened by
excluding randfold analysis.
# Load bash variables into memory
source .bashvars
# Append miRDeep2 to system PATH and set PERL5LIB
export PATH=$PATH:/home/shared/mirdeep2/bin
export PERL5LIB=$PERL5LIB:/home/shared/mirdeep2/lib/perl5
time \
${mirdeep2} \
${output_dir_top}/${collapsed_reads_mirdeep2_over17bp} \
${genome_fasta_dir}/${genome_fasta_no_spaces} \
${output_dir_top}/${mirdeep2_mapping_file} \
none \
${deep_dive_data_dir}/${mirbase_mature_fasta_no_spaces} \
none \
-t S.purpuratus \
-P \
-v \
-g -1 \
2>${output_dir_top}/miRDeep2-S.purpuratus-report.log# Load bash variables into memory
source .bashvars
tail -n 6 ${output_dir_top}/miRDeep2-S.purpuratus-report.logmiRDeep runtime:
started: 15:10:16
ended: 20:38:19
total:29h:28m:3sMiRDeep2 outputs all files to the current working directly with no way to redirect so want to move to intended output directory.
Output files will be in the format of result_* and error_*
# Load bash variables into memory
source .bashvars
for file in result_* error_*
do
mv "${file}" "${output_dir_top}/"
done# Load bash variables into memory
source .bashvars
ls -l | grep "^d"drwxr-xr-x 3 shedurkin labmembers 4096 Nov 2 08:39 03-Peve-lncRNA-dist_files
drwxr-xr-x 4 shedurkin labmembers 4096 Apr 22 13:18 10.1-Peve-sRNAsea-BLASTn-31bp-fastp-merged-cnidarian_miRBase_cache
drwxr-xr-x 4 shedurkin labmembers 4096 Apr 15 09:11 10.2-Peve-sRNAseq-BLASTn-31bp-fastp-merged_cache# Load bash variables into memory
source .bashvars
rsync -aP . --include='dir***' --exclude='*' --quiet "${output_dir_top}/"
rsync -aP . --include='map***' --exclude='*' --quiet "${output_dir_top}/"
rsync -aP . --exclude='mirgene*' --include='mir***' --exclude='*' --quiet "${output_dir_top}/"
rsync -aP . --include='pdfs***' --exclude='*' --quiet "${output_dir_top}/"
echo ""
echo "Check new location:"
echo ""
ls -l "${output_dir_top}/" | grep "^d"
FYI - eval=FALSE is set because the following command will only work
once…
ls --directory dir_* mapper_logs mirdeep_runs mirna_results* pdfs_*# Load bash variables into memory
source .bashvars
rm -rf dir_* mapper_logs mirdeep_runs mirna_results* pdfs_*
ls --directory dir_* mapper_logs mirdeep_runs mirna_results* pdfs_*ls: cannot access 'dir_*': No such file or directory
ls: cannot access 'mapper_logs': No such file or directory
ls: cannot access 'mirdeep_runs': No such file or directory
ls: cannot access 'mirna_results*': No such file or directory
ls: cannot access 'pdfs_*': No such file or directoryThe formatting of this CSV is terrible. It has a 3-column table on top of a 17-column table. This makes parsing a bit of a pain in its raw format.
# Load bash variables into memory
source .bashvars
head -n 30 "${output_dir_top}/result_22_04_2024_t_15_10_16.csv"miRDeep2 score estimated signal-to-noise excision gearing
10 6.4 1
9 6.2 1
8 6.1 1
7 5.7 1
6 6.2 1
5 9.1 1
4 7.2 1
3 4.6 1
2 3 1
1 2.3 1
0 1.9 1
-1 1.5 1
-2 1.1 1
-3 0.8 1
-4 0.7 1
-5 0.7 1
-6 0.7 1
-7 0.7 1
-8 0.8 1
-9 0.8 1
-10 0.8 1
novel miRNAs predicted by miRDeep2
provisional id miRDeep2 score estimated probability that the miRNA candidate is a true positive rfam alert total read count mature read count loop read count star read count significant randfold p-value miRBase miRNA example miRBase miRNA with the same seed UCSC browser NCBI blastn consensus mature sequence consensus star sequence consensus precursor sequence precursor coordinate
Porites_evermani_scaffold_1503_530575 110319.9 - 216386 216340 0 46 yes - - - - ucaggucuaggcugguuaguuu cuauaccagucuaguccu cuauaccagucuaguccuggcauguuucuugucaggucuaggcugguuaguuu Porites_evermani_scaffold_1503:46868..46921:+
Porites_evermani_scaffold_26_42156 46827.4 - 91841 89274 0 2567 yes - gga-miR-1467-5p_MIMAT0007345_Gallus_gallus_miR-1467-5p - - ucucagcucaccaaucucugcu cagggacuggugagcugauguc cagggacuggugagcugaugucauuuacugaucucagcucaccaaucucugcu Porites_evermani_scaffold_26:382571..382624:-
Porites_evermani_scaffold_910_418426 43145.7 - 84619 84496 0 123 yes - hsa-miR-33a-3p_MIMAT0004506_Homo_sapiens_miR-33a-3p - - caauguuucggcuuguucccg ggaacaagccgaaacauuuga caauguuucggcuuguucccguuuucgggaacaagccgaaacauuuga Porites_evermani_scaffold_910:118741..118789:+This will match the line beginning with provisional id and print to
the end of the file (represented by the $p. $ = end, p = print)
# Load bash variables into memory
source .bashvars
sed --quiet '/provisional id/,$p' "${output_dir_top}/result_22_04_2024_t_15_10_16.csv" \
> "${output_dir_top}/parsable-result_22_04_2024_t_15_10_16.csv"
head "${output_dir_top}/parsable-result_22_04_2024_t_15_10_16.csv"provisional id miRDeep2 score estimated probability that the miRNA candidate is a true positive rfam alert total read count mature read count loop read count star read count significant randfold p-value miRBase miRNA example miRBase miRNA with the same seed UCSC browser NCBI blastn consensus mature sequence consensus star sequence consensus precursor sequence precursor coordinate
Porites_evermani_scaffold_1503_530575 110319.9 - 216386 216340 0 46 yes - - - - ucaggucuaggcugguuaguuu cuauaccagucuaguccu cuauaccagucuaguccuggcauguuucuugucaggucuaggcugguuaguuu Porites_evermani_scaffold_1503:46868..46921:+
Porites_evermani_scaffold_26_42156 46827.4 - 91841 89274 0 2567 yes - gga-miR-1467-5p_MIMAT0007345_Gallus_gallus_miR-1467-5p - - ucucagcucaccaaucucugcu cagggacuggugagcugauguc cagggacuggugagcugaugucauuuacugaucucagcucaccaaucucugcu Porites_evermani_scaffold_26:382571..382624:-
Porites_evermani_scaffold_910_418426 43145.7 - 84619 84496 0 123 yes - hsa-miR-33a-3p_MIMAT0004506_Homo_sapiens_miR-33a-3p - - caauguuucggcuuguucccg ggaacaagccgaaacauuuga caauguuucggcuuguucccguuuucgggaacaagccgaaacauuuga Porites_evermani_scaffold_910:118741..118789:+
Porites_evermani_scaffold_72_87796 31859.5 - 62494 58655 1661 2178 no - spi-mir-temp-5_Stylophora_pistillata_Liew_et_al._2014_NA - - gagguccggacgguugaggguuauc caccccucauccaccaacuugaccucucu gagguccggacgguugaggguuaucaauuuauacuagucugcucaacuggaauuucugaaccaccccucauccaccaacuugaccucucu Porites_evermani_scaffold_72:198220..198310:+
Porites_evermani_scaffold_1503_530579 27800.7 - 54522 54405 0 117 yes - mmu-miR-710_MIMAT0003500_Mus_musculus_miR-710 - - ucaagucuaggcugguuaguuu cuacaccaguguagucuuggca cuacaccaguguagucuuggcaugcuucuugucaagucuaggcugguuaguuu Porites_evermani_scaffold_1503:47579..47632:+
Porites_evermani_scaffold_469_286650 23872.2 - 46837 46752 12 73 no - mco-miR-12071-5p_MIMAT0046816_Mesocestoides_corti_miR-12071-5p - - uaacgguuuguuucuuccacaaug ucuguccacuuugcgagccuugcugcug ucuguccacuuugcgagccuugcugcugaucgagccaacuuuugucuugcauaacgguuuguuucuuccacaaug Porites_evermani_scaffold_469:45430..45505:-
Porites_evermani_scaffold_594_331559 22726.6 - 44568 41387 0 3181 yes - nve-miR-2023-3p_MIMAT0009756_Nematostella_vectensis_miR-2023-3p - - aaagaaguacaagugguaggg cugccacuuguaucuucuuuca cugccacuuguaucuucuuucacguuuaucgaugaaagaaguacaagugguaggg Porites_evermani_scaffold_594:158195..158250:+
Porites_evermani_scaffold_253_199635 22568.8 - 44259 38319 0 5940 yes - osa-miR5539a_MIMAT0022175_Oryza_sativa_miR5539a - - uagaaaacucguguacgugacccu cgcacguagacgaguuuuuaac cgcacguagacgaguuuuuaacuuuugaugguuagaaaacucguguacgugacccu Porites_evermani_scaffold_253:200964..201020:+
Porites_evermani_scaffold_577_326309 20195.8 - 39607 39544 0 63 yes - bmo-miR-3329_MIMAT0015517_Bombyx_mori_miR-3329 - - ucauacaauggacaaggaucagg ugauucuuguucuuguuauaaua ugauucuuguucuuguuauaauaauuuuugucauacaauggacaaggaucagg Porites_evermani_scaffold_577:78358..78411:+This chunk provides a more convise overview of the data and it’s columns.
mirdeep_result.df <- read.csv("../output/11.1-Peve-sRNAseq-miRdeep2-31bp-fastp-merged-cnidarian_miRBase/parsable-result_22_04_2024_t_15_10_16.csv",
header = TRUE,
sep = "\t")
str(mirdeep_result.df)'data.frame': 2616 obs. of 17 variables:
$ provisional.id : chr "Porites_evermani_scaffold_1503_530575" "Porites_evermani_scaffold_26_42156" "Porites_evermani_scaffold_910_418426" "Porites_evermani_scaffold_72_87796" ...
$ miRDeep2.score : num 110320 46827 43146 31860 27801 ...
$ estimated.probability.that.the.miRNA.candidate.is.a.true.positive: logi NA NA NA NA NA NA ...
$ rfam.alert : chr "-" "-" "-" "-" ...
$ total.read.count : int 216386 91841 84619 62494 54522 46837 44568 44259 39607 35129 ...
$ mature.read.count : int 216340 89274 84496 58655 54405 46752 41387 38319 39544 34380 ...
$ loop.read.count : int 0 0 0 1661 0 12 0 0 0 0 ...
$ star.read.count : int 46 2567 123 2178 117 73 3181 5940 63 749 ...
$ significant.randfold.p.value : chr "yes" "yes" "yes" "no" ...
$ miRBase.miRNA : chr "-" "-" "-" "-" ...
$ example.miRBase.miRNA.with.the.same.seed : chr "-" "gga-miR-1467-5p_MIMAT0007345_Gallus_gallus_miR-1467-5p" "hsa-miR-33a-3p_MIMAT0004506_Homo_sapiens_miR-33a-3p" "spi-mir-temp-5_Stylophora_pistillata_Liew_et_al._2014_NA" ...
$ UCSC.browser : chr "-" "-" "-" "-" ...
$ NCBI.blastn : chr "-" "-" "-" "-" ...
$ consensus.mature.sequence : chr "ucaggucuaggcugguuaguuu" "ucucagcucaccaaucucugcu" "caauguuucggcuuguucccg" "gagguccggacgguugaggguuauc" ...
$ consensus.star.sequence : chr "cuauaccagucuaguccu" "cagggacuggugagcugauguc" "ggaacaagccgaaacauuuga" "caccccucauccaccaacuugaccucucu" ...
$ consensus.precursor.sequence : chr "cuauaccagucuaguccuggcauguuucuugucaggucuaggcugguuaguuu" "cagggacuggugagcugaugucauuuacugaucucagcucaccaaucucugcu" "caauguuucggcuuguucccguuuucgggaacaagccgaaacauuuga" "gagguccggacgguugaggguuaucaauuuauacuagucugcucaacuggaauuucugaaccaccccucauccaccaacuugaccucucu" ...
$ precursor.coordinate : chr "Porites_evermani_scaffold_1503:46868..46921:+" "Porites_evermani_scaffold_26:382571..382624:-" "Porites_evermani_scaffold_910:118741..118789:+" "Porites_evermani_scaffold_72:198220..198310:+" ...This provides some rudimentary numbers for the miRDeep2 output.
Further analysis is possibly desired to evaluate score thresholds, miRNA families, etc.
# Load bash variables into memory
source .bashvars
# Total predicted miRNAS
total_miRNAs=$(awk 'NR > 1' ${output_dir_top}/parsable-result_22_04_2024_t_15_10_16.csv \
| wc -l
)
echo "Total of predicted miRNAs: ${total_miRNAs}"
echo ""
# Matches to known mature miRNAs
mature_miRNAs=$(awk -F'\t' '$11 != "-" && $11 != "" {print $11}' ${output_dir_top}/parsable-result_22_04_2024_t_15_10_16.csv \
| wc -l
)
echo "Number of seed matches to known miRNAS: ${mature_miRNAs}"
echo ""
# Novel miRNAs
novel_miRNAs=$(awk -F "\t" '$11 == "-" || $11 == "" {print $11}' ${output_dir_top}/parsable-result_22_04_2024_t_15_10_16.csv \
| awk 'NR > 1' \
| wc -l
)
echo "Number of novel miRNAs: ${novel_miRNAs}"
Total of predicted miRNAs: 2616
Number of seed matches to known miRNAS: 2300
Number of novel miRNAs: 316