Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ENH: allow passing MAGs to the compute action #7

Merged
merged 3 commits into from
Aug 17, 2023
Merged

ENH: allow passing MAGs to the compute action #7

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
f9e9e50
ENH: allow passing MAGs to the compute action
misialq Aug 8, 2023
3b6de0c
flake8
misialq Aug 8, 2023
739b269
Update README
misialq Aug 8, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 7 additions & 0 deletions README.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,13 @@ You need to have QIIME 2 version 2018.4 or later. Also, regardless of which way

`conda install -c bioconda sourmash`

You will also need to install q2-types-genomics (unless your environment already has it):

```
conda install -c conda-forge -c bioconda -c https://packages.qiime2.org/qiime2/2023.5/tested -c defaults \
q2-types-genomics
```

To install the plugin, run the following command:

```
Expand Down
33 changes: 27 additions & 6 deletions q2_sourmash/_compare.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,28 +5,49 @@
#
# The full license is in the file LICENSE, distributed with this software.
# ----------------------------------------------------------------------------

from q2_types.distance_matrix import DistanceMatrix
from q2_sourmash._format import MinHashSigJsonDirFormat
import os
import re
import subprocess

import numpy
import skbio

def compare(min_hash_signature:MinHashSigJsonDirFormat, ksize: int, ignore_abundance: bool=True) -> skbio.DistanceMatrix:
from q2_sourmash._format import MinHashSigJsonDirFormat


def compare(
min_hash_signature: MinHashSigJsonDirFormat,
ksize: int,
ignore_abundance: bool = True
) -> skbio.DistanceMatrix:

np_file = 'tmp'
label_file = 'tmp.labels.txt'
command = ['sourmash', 'compare', str(min_hash_signature) + "/*", '--ksize', str(ksize), '-o', 'tmp']

command = [
'sourmash', 'compare', f'{min_hash_signature}/*',
'--ksize', str(ksize), '-o', 'tmp'
]
if ignore_abundance:
command.append('--ignore-abundance')

subprocess.run(' '.join(command), check=True, shell=True)

# load np_file as np.ndarray -> np_sim
np_sim = numpy.load(np_file)

# convert similarity to distance
np_dis = 1 - np_sim

# read labels into a list -> labels
labels = [os.path.basename(filename).strip().strip('.fastq.gz') for filename in open(label_file)]
labels = [
re.sub(
r'\.(fastq\.gz|fasta)$', '', os.path.basename(filename).strip()
) for filename in open(label_file)
]

# clean up
os.remove(np_file)
os.remove(label_file)

return skbio.DistanceMatrix(np_dis, labels)
77 changes: 59 additions & 18 deletions q2_sourmash/_compute.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,40 +5,81 @@
#
# The full license is in the file LICENSE, distributed with this software.
# ----------------------------------------------------------------------------

from q2_types.per_sample_sequences import SingleLanePerSampleSingleEndFastqDirFmt, FastqGzFormat
import qiime2.util
import pandas as pd
from q2_sourmash._format import MinHashSigJsonDirFormat
import glob
import os
import subprocess
import glob
import sys
from typing import Union

import pandas as pd
import qiime2.util
from q2_types.per_sample_sequences import (
SingleLanePerSampleSingleEndFastqDirFmt
)
from q2_types_genomics.per_sample_data import MultiMAGSequencesDirFmt

from q2_sourmash._format import MinHashSigJsonDirFormat


def compute(sequence_file:SingleLanePerSampleSingleEndFastqDirFmt, ksizes: int, scaled: int, track_abundance: bool=True) -> MinHashSigJsonDirFormat:
def _duplicate_mag_seqs(manifest, sequence_file):
output = MinHashSigJsonDirFormat()
for i, row in manifest.iterrows():
_, bin_id, src_fp = row["sample-id"], row["mag-id"], row["filename"]
src_fp = os.path.join(str(sequence_file), src_fp)
dest_fp = os.path.join(str(output), f'{bin_id}.fasta')
qiime2.util.duplicate(src_fp, dest_fp)
return output

#read in FastqManifestFormat to convert from sample name to filename
manifest = pd.read_csv(os.path.join(str(sequence_file), sequence_file.manifest.pathspec), header=0, comment='#')

def _duplicate_seqs(manifest, sequence_file):
output = MinHashSigJsonDirFormat()
for seq_file in glob.glob(os.path.join(str(sequence_file), '*.fastq.gz')):
src_fp = str(seq_file)
filename = os.path.basename(src_fp)
sample_id = list(
manifest[manifest['filename'] == filename]['sample-id']
)
dest_fp = os.path.join(str(output), f'{sample_id[0]}.fastq.gz')
qiime2.util.duplicate(src_fp, dest_fp)
return output


def compute(
sequence_file: Union[
SingleLanePerSampleSingleEndFastqDirFmt, MultiMAGSequencesDirFmt
],
ksizes: int,
scaled: int,
track_abundance: bool = True
) -> MinHashSigJsonDirFormat:
# read in FastqManifestFormat to convert from sample name to filename
manifest = pd.read_csv(
os.path.join(str(sequence_file), sequence_file.manifest.pathspec),
header=0, comment='#'
)

for seq_file in glob.glob(os.path.join(str(sequence_file), '*fastq.gz')):
filepath = str(seq_file)
filename = os.path.basename(filepath)
sampleid = list(manifest[manifest['filename']==filename]['sample-id'])
# print(sampleid, filename)
qiime2.util.duplicate(filepath, os.path.join(str(output), sampleid[0]+'.fastq.gz'))
if isinstance(sequence_file, MultiMAGSequencesDirFmt):
output = _duplicate_mag_seqs(manifest, sequence_file)
else:
output = _duplicate_seqs(manifest, sequence_file)

command = ['sourmash', 'compute', str(output) + "/*", '--ksizes', str(ksizes), '--scaled', str(scaled)]
command = [
'sourmash', 'compute', f'{output}/*',
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

in our next major release of sourmash, sourmash v5, we will be removing sourmash compute per sourmash-bio/sourmash#1286.

It's been replaced with sourmash sketch; the following should work:

ksizestr = [ f"k={k}" for k in ksizes ]
ksizestr = ",".join(ksizestr)
abundstr = "noabund"
if track_abundance:
   abundstr = "abund"

command = ['sourmash', 'sketch', 'dna', f'{output}/*', '-p', f'scaled={scaled},{ksizestr},{abundstr}'

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good to know, thanks for info! Should we then update it in a separate PR, once v5 is out?

'--ksizes', str(ksizes),
'--scaled', str(scaled)
]

if track_abundance:
command.append('--track-abundance')

subprocess.run(' '.join(command), check=True, shell=True, cwd=str(output))

for seq_file in glob.glob(os.path.join(str(output), '*fastq.gz')):
os.remove(seq_file)
if isinstance(sequence_file, MultiMAGSequencesDirFmt):
ext = '*.fasta'
else:
ext = '*.fastq.gz'
for fp in glob.glob(os.path.join(str(output), ext)):
os.remove(fp)

sys.stdout.flush()

Expand Down
48 changes: 25 additions & 23 deletions q2_sourmash/plugin_setup.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,19 +5,15 @@
#
# The full license is in the file LICENSE, distributed with this software.
# ----------------------------------------------------------------------------


import q2_sourmash

import qiime2.plugin
from qiime2.plugin import Plugin, Metadata, Str, List, Citations, SemanticType, TextFileFormat, ValidationError
from qiime2.plugin import model
import qiime2.util
from q2_sourmash._compute import compute
from q2_sourmash._compare import compare
from q2_types.distance_matrix import DistanceMatrix
from q2_types.sample_data import SampleData
from q2_types.per_sample_sequences import SequencesWithQuality
from q2_types.sample_data import SampleData
from q2_types_genomics.per_sample_data import MAGs
from qiime2.plugin import Plugin, Citations

from q2_sourmash._compare import compare
from q2_sourmash._compute import compute
from ._format import MinHashSigJsonDirFormat, MinHashSigJson
from ._types import MinHashSig

Expand All @@ -29,8 +25,8 @@
citations=Citations.load('citations.bib', package='q2_sourmash'),
description=('This QIIME 2 plugin wraps sourmash and '
'supports the calculation and comparison of '
'minhash signatures.'),
short_description='Plugin for generation of minhash signatures.'
'MinHash signatures.'),
short_description='Plugin for generation of MinHash signatures.'
)


Expand All @@ -43,20 +39,26 @@

plugin.methods.register_function(
function=compute,
inputs={'sequence_file': SampleData[SequencesWithQuality]},
parameters={'ksizes': qiime2.plugin.Int,
inputs={'sequence_file': SampleData[SequencesWithQuality | MAGs]},
parameters={
'ksizes': qiime2.plugin.Int,
'scaled': qiime2.plugin.Int,
'track_abundance': qiime2.plugin.Bool},
'track_abundance': qiime2.plugin.Bool
},
outputs=[('min_hash_signature', MinHashSig)],
name = 'compute sourmash signature',
description = 'Computes a sourmash MinHash signature from fasta/q files.'
name='Compute sourmash signature',
description='Computes a sourmash MinHash signature from fasta/q files.'
)

plugin.methods.register_function(function=compare,
inputs={'min_hash_signature':MinHashSig},
parameters={'ksize': qiime2.plugin.Int,
'ignore_abundance': qiime2.plugin.Bool},
plugin.methods.register_function(
function=compare,
inputs={'min_hash_signature': MinHashSig},
parameters={
'ksize': qiime2.plugin.Int,
'ignore_abundance': qiime2.plugin.Bool
},
outputs=[('compare_output', DistanceMatrix)],
name = 'compare sourmash signatures',
description = 'Compares sourmash signatures and calculats Jacaard distance matrix.'
name='Compare sourmash signatures',
description='Compares sourmash signatures and '
'calculates Jacaard distance matrix.'
)