soedinglab / plass

sensitive and precise assembly of short sequencing reads
https://plass.mmseqs.com
GNU General Public License v3.0
132 stars 14 forks source link

Using already predicted ORFs from reads #24

Open genomewalker opened 4 years ago

genomewalker commented 4 years ago

Hi I would like to use the ORFs that I predicted using another approach. I checked the documentation and I haven't seen an option to deactivate the extractorfs step. Is this possible or do you have any suggestions on how to use already predicted ORFs?

Thank you very much Antonio

martin-steinegger commented 4 years ago

For now there is only a hacky ways to todo this. You could start your job and immediately stop it. Now call:

plass createdb yourorfs.faa orfs 
touch tmp/latest/nucl_6f_start 
touch tmp/latest/aa_6f_start 
touch tmp/latest/nucl_6f_long
touch tmp/latest/aa_6f_long
ln -s orfs tmp/latest/aa_6f_start_long 
ln -s orfs.index tmp/latest/aa_6f_start_long.index
ln -s orfs.dbtype tmp/latest/aa_6f_start_long.dbtype
ln -s orfs_h tmp/latest/aa_6f_start_long_h
ln -s orfs_h.index tmp/latest/aa_6f_start_long_h.index
ln -s orfs_h.dbtype tmp/latest/aa_6f_start_long_h.dbtype

Now restart your run by calling the initial plass command.

genomewalker commented 4 years ago

Thanks Martin!

I have a couple of hundreds of metagenomes to assemble. For now, I will try to modify the PLASS assembler workflow and add these steps and see if it works.

Antonio

genomewalker commented 4 years ago

Hi Martin I modified assembledb.sh to include the hack you suggested and seems to works correctly. I can reduce the number of assembled proteins with the default NN filter value.

if notExists "${TMP_PATH}/nucl_6f_start"; then
    touch "${TMP_PATH}"/nucl_6f_start
fi

if notExists "${TMP_PATH}/aa_6f_start"; then
    touch "${TMP_PATH}"/aa_6f_start
fi

if notExists "${TMP_PATH}/nucl_6f_long"; then
    touch "${TMP_PATH}"/nucl_6f_long
fi

if notExists "${TMP_PATH}/aa_6f_start_long"; then
    touch "${TMP_PATH}"/aa_6f_long
    ln -s "${INPUT}" "${TMP_PATH}"/aa_6f_start_long
    ln -s "${INPUT}".index "${TMP_PATH}"/aa_6f_start_long.index
    ln -s "${INPUT}".dbtype "${TMP_PATH}"/aa_6f_start_long.dbtype
fi

if notExists "${TMP_PATH}/aa_6f_start_long_h"; then
    ln -s "${INPUT}"_h "${TMP_PATH}"/aa_6f_start_long_h
    ln -s "${INPUT}"_h.index "${TMP_PATH}"/aa_6f_start_long_h.index
    ln -s "${INPUT}"_h.dbtype "${TMP_PATH}"/aa_6f_start_long_h.dbtype
fi

Some numbers:

file         format  type     num_seqs  sum_len  min_len  avg_len  max_len
orfs_aa.faa  FASTA   Protein     2,912   82,008       17     28.2       84
orfs_nt.faa  FASTA   Protein     4,236  130,973       21     30.9       91
raw.faa      FASTA   Protein     5,609  173,304       21     30.9       92

Where orfs_aa are the predicted ORFs in aa, orfs_nt are the predicted ORFs in nt and raw are the raw reads.

These are ancient reads (short and damaged), and using the default workflow, I was getting many proteins that I believe might be spurious.

I will rerun all metagenomes with the modified workflow and compare the outcome in terms of annotations and other parameters I calculate. I keep you posted.

Thank you very much Antonio

martin-steinegger commented 4 years ago

Thank you for sharing the data. I think the spurious proteins might be introduced by our change to add complete short Orfs. Even if they were not assembled. E.g. if the is a complete Orf encoded on one read then we add them to the final result. We try to filter these out with the NN, but the NN seems to not always be able to remove them. Do you see that many of the spurious proteins are short?