Given that you have a very low fraction of bacterial reads, which is a common problem in the field, I think a useful contribution from your data would be to create a panel of primers to amplify community members that you see are present. This would give you more resolution than 16s but allow you to avoid more of the host sequencing data. However, the usefulness of such a panel would be bounded by how it would be adopted by others in the field. It would probably be most useful if you applied it to this fish farm repeatedly, but I'm not sure if doing so is biologically interesting.

Since you have the resolution to go below phylum, I think it would be interesting to focus on that more in the discussion.

I agree, but I think this would be worth mentioning in the abstract, and perhaps in the last paragraph of the introduction, to better prepare your readers about the types of results you are going to present

The percentage or fraction? If percentage, less than 1% is incredibly small and I would question any results in this report. How many reads total was this? If you used a very high depth, you might capture a substantial portion of the community.

I imagine you might have a huge amount of drop out here by applying kraken first and then assembling with megahit. I would either: 1. Map reads to carp first, and then assemble anything that doesn't map 2. Assemble everything and then filter out carp contigs.

Why not the species level?

It might be worth mapping back to the host genome if you have one prior to performing taxonomic classification.

I would also be interested to see nonpareil curves of your sequencing data before and after host mapping. I would be curious if you reached saturation of the community -- this can usually be better assessed with raw sequencing data than with taxonomically classified reads.

Can you add what functions you used to do this in the Vegan package?

What filters did you use here? I'm curious how many reads were lost to filtering.

Do you have any idea if the bacterial load of the water, and therefore the skin, of the carp was much higher than for fish observed in the wild, or typically sequenced with 16s? I'm wondering if there was more bacteria than usual, and that was why you were able to get enough bacterial reads to perform an analysis

Would you be able to provide some specific examples of the economic importance (even half a sentence)? I'm not a carp expert so I have no idea what these might be!

Can the microbiome also be pathogenic for carp?

I think this section doesn't highlight the massive challenge in trying to get shotgun metagenomic sequencing data from fish. In the experiments where we have tried (killifish, different tissues), we end up with 98 or 99% killifish (host) reads. 16s allows us to amplify and get just the microbial signal.

We have talked about trying to do a more balanced marker gene panel, but that has methodological problems like not having as many tools and determining the best marker genes to use.

It would be nice if these challenges were better represented. I think the reason this gap exists is methodological (hard to get shotgun sequencing from fish), not for lack of interest

I'm curious the extent to which these studies investigated farm vs. wild fish, and if you think that would make a difference on microbiome. It might be helpful to include that distinction in when covering literature in the introduction, given that you see some results that you don't expect relative to other observations in the field.

I may have missed it, but I didn't see a methods section for peptide discovery/annotation. Would you be willing to add this?

Is it possible to be more specific about the sequencing data that would be needed to answer this question? I think that could be value-added for the community to have this clearly spelled out

would these be post-assembly errors or assembly errors? If post-assembly, can you add details about what this means?

It isn't clear in this section whether you ran BUSCO on the genomes in genome or protein mode. I think it would be helpful to see both sets of results -- the genome alone would tell us how many single copy genes are in the assembly, while the proteins compared with this will tell us how well annotation went for these genomes as well.

Are your methods versatile enough to allow this approach for new species pairs if there is enough RNAseq data? What is the minimum number of RNA seq samples needed?

Is there a way to look at this more systematically? This is hugely valuable.

Do you think this is because of the expression levels of genes or some other signal? I think unpacking this more would be a big value-added for understanding this model.

This is very cool. it would be interesting to correlate this with the number of failed clinical trials for therapies developed in mouse and applied to humans. It might also be interesting to see if there are other DNA signals that could be used to improve DNA performance. Promoter sequences come to mind, but there might be others.

Then what was it due to?

What were the other, nearer embedding neighbors? copy number variants of the same genes? genes in the same pathway?

I don't understand the importance of this note as currently written.

Do you have confirmation that all information capture in the embeddings is synonymous between these 1:1 orthologs? For example, doe the orthologs share the same promoters and transcription factor still? I think both of these would be fairly straightforward to check at scale.

Thank you for including these details, this is very helpful

I don't think this paragraph does justice to the work already undertaken and I don't think it highlights why there was a gap for this work to fill. I think the goals of the previous compendia that you included were very different than the goal of this paper, and I think that's ok! Unpacking that a bit more in this intro would be useful. Otherwise, it makes it sound like previous researchers made mistakes and that's why we need this paper, which I don't think is true.

Old URL I think (though it still re-directs to the correct place!)

Thanks for putting this together! I took a look at the repository and noticed a few changes that I think could help make CELEBRIMBOR more user friendly. 1. Would you be willing to add tool versions to your environment.yml file? This will help make sure this will still be installable over time, and help the docker container match in results to users who deploy this on e.g. an hpc 2. Does create_plots.py belong in the scripts/ directory instead? 3. The hashes after the second equals sign in some of your yaml files will make the environments difficult to install across different operating systems (linux vs. mac) (ex https://github.com/bacpop/CELEBRIMBOR/blob/main/envs/Snakemake.yaml) 4. The readme uses the old repo name in the clone/cd instructions 5. Lastly, have you explored using something like a click interface & making the tool conda-installable? I know this is a big lift, but in my experience it makes it so much easier for others to pick up and use, including dropping the pipeline into larger pipelines. This pipeline might provide some inspiration for how to accomplish this: https://github.com/metagenome-atlas/atlas

can you be more specific here? which functions in these packages accomplish the task (e.x. I believe it is vst() for DESeq2). Can you also cite DESeq2 instead of DESeq? I think DESeq has been retired.

would you consider switching this plot to an upset plot (R packages upsetr or complexupset) instead of a venn diagram? For many intersections, upset plots are a bit easier to understand than venn diagrams

are these number right? are 4 samples from one subject? 113*2=226, not 228

Did you do anything to account for these things in your analysis?

I'm so glad you included this! I have struggled with this a lot in my own research so I'm so glad to see it explicitly mentioned here.

This is very nice.

Is this supposed to be a less than sign instead of a greater than sign?

can you prepend this with "cancer-type-specific" so that its clear inline what this means without having to prematurely jump to a future section?

Can you provide a little more context as to what this means? Are all samples consistently analyzed or is there some normalization that takes place as well?

I get a "Page not found" error when i navigate to this URL. Is the repo still private or is there a typo in the URL? I would love to look at/give feedback on the code as well!

This sentence switches to passive from active voice, and the next sentence is active again ("We"). Would it be possible to make it active voice as well? without that, it doesn't sound like you ran the GSEA analysis and got it from somewhere else, which is a little confusing

Meaning short read Illumina bulk RNA-seq data?

Is this the recommended amount of RAM to run chlomito? this is quite high.

This link gave me a 404

I'm assuming the second metric relies on mapped reads. Did you consider identifying spanning reads as further evidence for your tool? If a read spans an organellar genome sequence and nuclear genome sequence (perhaps with k=21 bp overlap at minimum, or potentially higher), then I think that would show evidence of an HGT event

Is chlomito available on GitHub? Where can the docker image be downloaded from?

Would you be willing to provide details on the quality of these two genomes? How well known are the chloroplast and mitochondrial sequences in these (do they have gold-standard labels?)?

Can you provide which date you accessed this data?

I'm curious about the method employed here and the rationale behind it. 1. What was the size of the sliding window? So that mean nucleotide at position X in a genome would be captured many times, each time multiple times for a sliding window, capturing all 3k bases before and after it in a series of different windows? 2. Why eliminate sequences less than 200bp in length if they are "real"? Many non coding RNAs or peptides are encoded by short sequences. Does this decision limit your embedding space to sequences > 200 bp long?

It would be super helpful at this point to breakdown a bit more they type of data this model was trained on -- genomes, transcriptomes, thinks in GenBank, etc.

Would it make sense to include citations to Nucleotide Transformer and DNABERT as well? Unless I'm misunderstanding the application space, these seem like seminal efforts in this area.

I'll post this on github as well, but are these models available for use? I know I could use UniDL4BioPep and their models, but I think your tool and would love to be able to use the models you built!

It would be interesting to unpack this sentence a bit more in the discussion. What do you mean by this, and can you provide any references for this?

Why do you think performance is so bad for blood brain barrier?

The way this is plotted is really confusing. I would suggest making the handcrafted models as a vertical line or a dot instead of its own bar. The reference line approach would make it clear that you're not re-running those models in some way like you are the other two.

This is sort of a dumb question at this point, but are all classifiers explicitly binary? How would you recommend a user who is interested in many classes of peptide bioactivity (or who plans to say, predict all peptides in a metagenome and then predict bioactivities) go about assessing all classes of bioactivity?

This is nice :)

Would it be possible for a user to provide some negative peptides, but then to also have those negative peptides supplemented by those generated in Step 2? This might be nice for users that have small input data sets if it's possible.

It would be super helpful if you could provide estimates for the minimum number of required sequences or any other information you can add to help a user build intuition for how many sequences they need to provide

I'm curious if using degenerate amino acid alphabets (dayhoff encoding, hydrophobic-polar, etc) would further improve classification accuracy or show interesting patterns.

it would be nice to see how many peptides were dropped at each stage of filtering. I'm also curious since you dropped so much data if the model would be more generalizable if that data were included somehow.

would you be willing to add a license to the repository so terms of re-use of your work are clear?

Sorry if I missed this, but can you report the total size of each group?

typo i think :)

curious if you tried without these groupings -- ie, how dissimilar are some of the peptides that were placed into combined groups, and could the model have done well without these groupings

Similarly, how often does this occur?

How often does this occur?

Did you try this as a multi-classification problem and end up with a better performance with binary classifiers? Or was the underlying model you used limited to binary classification?

I'm curious if this task is different than determining whether an amino acid sequence of 2-50 aas is a bioactive peptide, or if one must first know that the sequence is a peptide to then apply these tools to categorize the peptide sequence into a functional class.

One question that this name, and the abstract in general, left me with is whether this method is extensible beyond microbial fermentation peptides. I will continue to read to find out, but I'm wondering if another sentence might be added to clarify this.

I'm also curious if microbial fermentation peptides include all known classes of peptides, or if there are some functional classifications that might not be labelled by CICERON because CICERON has not seen them before. Again I will continue reading to hopefully find out!

have you started work around this? have you had any community buy in? do you have realistic timelines and goals around achieving this?

This is a massive drawback and essentially makes the addition of all of this provenance info moot at worst and irrelevant to many FASTA files at best.

Would it be possible to design a metadata standard that didn't rely on a header that would make most tools unable to process the data? Could it be extensible to other types of FASTAs?

double parentheses typo

Can you expand on what you mean by this?

missing comma

Will new tools be compatible with these legacy features? I'll be curious as I continue reading the paper whether you have tried using a FASTA with this header with popular tools (BWA, seqkit, seqtk, samtools, etc)

Do you need buy in from each of these communities for your metadata standard to be a success? How do you plan to get that buy in?

Does NCBI fundamentally not allow for these metadata fields, or are they not inserted by the users who upload the data? I think creating a new metadata standard (as presented here) while in theory could solve some of these issues, compliance by those who upload genomes will always be an issue no matter what standard is used. I think researchers default to not including information when they are unsure about that information, or unsure of themselves at time of upload, which has historically been a rather stressful process.

typo :)

These are really great examples!

I think NCBI handles much of this by putting the information directly in the FASTA header for each contig. The accession itself given to each contig creates a link but I think it does act as a small tracker of provenance

These citations follow a different format

Would you be willing to use a synonym here? After reading the abstract and first paragraph, I'm searching for clues that this metadata standard might also apply to other types of sequencing data (other FASTAs with e.g. assembled transcripts, amino acid sequences, genes, etc), and seeing transcription here and not referring to the central dogma is a little distracting

how did you deal with shared KOs between pathways when doing pathway level analysis?

many of these sequences contain detectable human sequences. I would be curious for you to run the human genome against your databases and see what functional profile is returned. it would let users know whether they need to do host filtering before applying this approach. if they didn't need to take that step (or any other QC), that would be huge time savings

how did you come up with parameters, and how do you know they are the best to use? how would you advise users to choose between k-mer sizes for their own applications?

reiterating the point from above, how does your approach break down with increasing evolutionary divergence from the reference, and how is that different from other tools. Soil might be a good ecosystem to test drive this in, and I think the CAMISIM tool allows you to introduce mutations from a reference in a known ratio/identity etc

I did some work similar to this with the pfam database a couple years ago: https://github.com/taylorreiter/2021-pfam-shared-kmers

I'm curious if you did any sort of analysis to see if there is shared kmer content between orthologous groups, or if high shared content (as is observed in pfam) would limit the ability of this approach to be generalized to other databases.

shorter

can you do these same analyses with a tool like HUMANN2 or something else that is typically used to do functional profiling and compare the results? can you show that you capture more functional units than other tools, or is your method only faster? would you need additional database above just KEGG to make the comparison fair, and is that possible with the approach you have outlined here?

I think the 2019 HMP IBD paper has a supplemental figure where they have KOs for each sample. it would be interesting to compare against those results for those samples to see if you get the same or different results (super set, subset, etc).

is the code for this also in the linked github repo? I couldn't find it, but I think it's an interesting application. It would be nice if something similar could be implemented for sourmash taxonomy results

pairwise distances between KOs obtained by comparing sourmash sketches, right?

I think this is too strong of a statement for the results presented. What about divergence in the metagenome vs. what's in the database? while using an amino acid k-mer will overcome some of this, I would expect diamond to better capture functional potential of a metagenome when the genomes are not in reference databases (I haven't explicitly done this test though so I don't know).

again, I think this comparison is unfair since you aren't using assembled genomes

this is not always a benefit. The "alignments" output by diamond can be super useful if the user wants to go back and do a targeted alignment of a specific gene of interest.

what about wall time? because diamond can be threaded, which is a huge plus, while sourmash gather cannot

this doesn't make a lot of sense as an application though right? kofamscan is designed to run on ORFs predicted from assembled genomes, not on metagenome reads?

I noticed that this repo doesn't have any unit tests and that the python script only contains 58 lines of code. Would it be possible to include this approach directly in sourmash?

even more than this, many of these algorithms are limited to the setting of assembled (meta)genomes, and there are a substantial number of studies showing that short read assembly often fails for metagenomes, especially for those from complex communities. If your method can work directly on short reads, I think that is a huge strength that is worth highlighting

(I believe DIAMOND-based approaches will also work quite white on short reads, but many of the others do not. while I have used diamond metagenomes against small databases [see the serratus rdrp paper for inspiration here], I'm not sure how well it would scale to whole metagenomes against all of e.g. KEGG).

Aren't they called KEGG Orthologs, which is abbreviated to KOs?

This statement isn't clear to me. It seems like the KOALA and kofamSCAN algorithms do consider orthology, can you expand this statement to make it clear what this means?

can you provide citations for this point, both before and after the semi colon?

last or lowest common ancestor?

I have some questions about this statement. 1. do you know if there is enough non-human single cell data to do a similar study on a different organism? Perhaps mouse might have enough for example? 2. Do you imagine this approach could be used for cross-species analysis, for example in which a mouse is compared to a human, or do you think it is limited to within species analysis?

This is a very clever manipulation of the input data that allows it to be analyze in a new way. As someone who is relatively new to this field, it would be very helpful if either in this section or the introduction if you could provide references for any approach within sequencing data that does something similar. If there is no similar approach to date, that would also be helpful to highlight. I think the idea of using an embedding is fairly common (e.g. word2vec), but it would be helpful to know the boundaries of innovation for this particular approach.

I'm really excited by the approach profiled in this paper. I think it's a very clever use of chemical reactions and stoichiometry. However, I'm concerned about the pre-processing step mentioned here. Do you have a sense of how lossy SUPER-FOCUS is, especially for soil microbiomes? Typically, metagenome analysis of soils can lose up to 80% of information due to the system complexity and the amount we have yet to observe. Depending on the fraction of loss of functional information, how do you expect that to impact that results presented in this study?

Did you explore the option of merging the StringTie and Trinity transcriptomes? Tools like orthofuser or evidentialgene could be used for this. I would be curious how well this might work. Can you think of any reason this would be a bad idea?

It would be super helpful in this section if you could compare the transcriptome assembly sequence themselves -- similarity and containment. I think this would help build intuition for the difference in transcriptome content between the different methods. A tool like sourmash can be used to estimate these metrics (sourmash sketch, and then sourmash compare. Then you can use sourmash plot to visualize). I'm curious if the trinity assembly has lower mapping because it has fragmented transcripts which are resolved with the use of the genome.

small typo -- period should be a comma :)

Can this be turned off? I can imagine for e.g. metagenomics one might not want these removed.

Just curious if you've thought of benchmarking with T2T assemblies. I think it could be useful and cool (although perhaps beyond the scope of this preprint) to see how much contam changes based on the completeness of the contam reference.

I'm also curious if these genomes needed to be masked at all (repeats, or ribosomal or other highly conserved sequences) so that off-target mapping doesn't occur. https://www.seqanswers.com/forum/bioinformatics/bioinformatics-aa/37175-introducing-removehuman-human-contaminant-removal

Can you change some of the early language to reflect this? The bbduk method is very different from mapping -- especially depending on the k-mer size that you end up using. It also will have very different output files (no BAM!), which is great, as BAMs are huge and a lot of people like to avoid them.

This figure is a bit blurry, would it be possible to update to a higher resolution image? Additionally, this overview is much more useful than the one currently in the GitHub readme. Would you be willing/able to swap out the README figure and place this one in instead?

In this manuscript, the authors present shournal, a tool to help with tracing shell commands that have been run on Linux computers. shournal sits in a space between iterative computational experiment and codifying those steps in a workflow. I'm excited by the concept of radical repeatability that lightweight tools like shournal could usher in.

I was unable to install shournal from the instructions on the github page, so this review does not cover feedback on the tool itself. I was eager to try on the snakemake integration and was sad not to be able to. I tried to install on an AWS EC2 instance (Ubuntu, t2.micro, using the latest release of shournal).

From a high level usability and adoption perspective, I think two things currently decrease the likelihood of shournal's broad adoption. First, the fact that there is no mac or windows distribution decreases shournal's audience. Second, the fact that shournal may be ineffective on HPCs further limits the audience (both by shournalk and by the event history not tracing over multiple machines These limitations do not decrease it's conceptual addition to the field, but will decrease the likelihood of adoption.

Culturally, I think there are pros and cons to shournal. On the pro side, I think having more tools in the reproducibility arsenal is a positive thing. Shournal meets scientists where they're at as they determine the best scripts to run on their data. However, I worry that reliance on shournal could lead to sub-par documentation for computational experiments. If researchers are in the habit of recording their commands and with notes, reliance on shournal may change this process, removing helpful metadata from command recordings. It is difficult to know how a tool like shournal could change the overall working habits of e.g. bioinformaticians, but it would be interesting to conduct a study on how adoption of shournal improves or detracts from reproducibility and documentation. (to be clear I am definitely NOT suggesting that that be done as part of this paper! But I think shournal could encourage a seachange in computing documentation, so it would be interesting from a metascience perspective to understand the benefits and drawbacks of those changes, and then how shournal could eventually be modified to reduce the drawbacks.)

One of the limitations suggested in the supplement is that, "provenance of binary executables is not tracked." Would it be possible to parse the help messages of binary executables or look at the stdout for version numbers or other tells of the software? This is far more inelegant than shournalk's current approach, but I wanted to supply it as a brainstorming idea in case the authors find it useful to iterate from. Alternatively, could the checksum of the binary executable be tracked?

Lastly, I left comments inline on the manuscript itself, but I also wanted to note that the first paragraph of the supplement provides important background knowledge that I think would be better served in the introduction of the paper if there is space to include it.

I think this URL has a typo in it

Clever -- so something like ls or cd would be totally ignored, but any program that actually looks at the data of a file will register, rigth?

This is slightly confusing wording for me -- does this mean typical shournal runs? For me, workflow is conflated with workflow engines, which I typically assume have a fairly substantial overhead (e.g. a snakemake workflow with a DAG of over 1m processes can take ~16gb of RAM to run)

What does this mean? It seems like shournal still needs to be activated at each shell session. Are you recommending that this be achieved by adding activation to the bash profile/bash rc, or is this only a statement that shournal was engineered to have a low overhead such that it could run indefinitely, even simulataneously with processes that demand high CPU, i/o, and/or RAM?

Do you have space to expand on this concept, even with another half sentence or sentence, on how the goals of shournal differ from the previously mentioned tools? I'm a huge nextflow & snakemake user and a big proponent of repeatable and reproducible computation, but I don't have the background to contextualize this comment which I think draws away from the potential impact of shournal in this space.

I think this citation may be broken

by virtue of being recorded so a scientist can go back and re-trace their steps, or through some other mechanism?

Is the unit for the y axis seconds?

For clarification, this is only if the original scripts or configuration files have been deleted, not if the data files (like a FASTQ or something) have been removed?

Can you make the code you used in this preprint available here as well?

were you time-limited in examining these or are they more difficult to tell than the above cases in whether the results represent true biology?

It would be helpful to know whether: 1. All of the data and code were available to try to exactly repeat what Poore et al did. If not, what ingredients are missing? 2. If you are able exactly repeat what Poore et al ran, do you get the exact same results? If not, is it because they didn't report e.g. a random seed value, or does it seem like the code that is reported isn't what was actually used?

Given that human contamination is a leading issue here, I think it could be interesting to show what mapping against a human pangenome accomplishes in terms of reduced number of human sequences in unmapped reads. Similarly, a program like bbmap's bbduk.sh could be used to remove human sequences. It's similar to kraken in that it detects matches based on k-mers. While neither of these things are necessary to prove the points in this preprint, I think this preprint might receive a lot of attention and I think it could be a gift to the community to demonstrate the most effective ways to reduce human reads in a sample when that is the goal. This would also have applications to the metagenomics field for host-based samples, where removal of human reads from e.g. gut microbiomes should be performed before data deposition.

Thanks for the preprint! It's very cool that enough public data exists to do this sort of analysis. I have some concerns about the "Imputed Metagenomics" section. I think the success of the picrust approach is dependent on how well the species in the microbiome are represented in reference genomes. I have some doubts that corn microbiome species would be well represented in public databases, especially for ecosystem-specific genes in the pangenome. I think it would be a huge value-added and validation if you were able to show that picrust accurately predicts the function of corn metagenomes. I think you could do this either by identifying a publicly available paired data set for each of your sample types (soil, rhizosphere, roots, and leaves) and doing functional annotation on the shotgun metagenome versus picrust on the 16s. If this doesn't exist, you could use a tool like groupM to pull out 16s sequences from a metagenome and analyze those with picrust, and then compare it to the functional annotation of the metagenome. Without this type of analysis, I think the "Imputed Metagenomics" section is overstated.

Thanks for the lovely paper! I have some comments on the section, "Binning of Prokaryotic Genomes and Bin Refinement." First, I'm curious if you tried to decontaminate any bins. I think this could be an interesting step, whether they initially had contamination greater than or less than 5%. I think GUNC can be used for this. I've also contributed to a software called charcoal that could work for this. Second, I have concerns about the bin refinement you performed. I think assemblers (esp short read assemblers, so may be less relevant) typically break contiguous sequences and produce fragments when there is either incomplete sequencing or strain variation that causes a bifurcation in the underlying graph. Could your refinement technique be merging contigs together that would never be combined in nature? It could be interesting to investigate the underlying assembly graph for some of these merges, either in the short or long reads. I see you think about this in the paragraph, "This bin refinement method allows the merging of contigs..." I think it might be worth additional investigation to ensure this is not happening, as the methods you use here could be precedent-setting for long read metagenomes (esp those sequenced without accompanying hi-c sequencing data).

Thanks for this effort, I really enjoyed your preprint. I loved the note at the beginning about the length of the preprint :)

I'm very excited by mTMAT. In the past, I've struggled to find robust statistical techniques to analyze longitudinal microbiome data, so I'm very excited to have another method in my toolkit. After reading your preprint, I have a few questions about the versatility of mTMAT. 1. Would it be possible to use mTMAT on shotgun metagenomic data? I'm envisioning longitudinal samples for which species or strain counts have been inferred. 2. How would mTMAT behave in the face of large microbial transitions over time. For example, what about microbiomes where the composition fluctuates diurnally? Or, a large change happens over the course of a year (like the onset of inflammatory bowel disease) 3. The language at the beginning of the paper emphasis association with disease state. I'm curious if mTMAT is flexible for application to any variable set. For example, could I apply this tool to wine microbiomes? I assume yes from the way the pregnancy data set is presented but wanted to check.

zDB looks really cool, I'm excited to try it! Have you thought about including an upset plot instead of a venn diagram? I think it could help make intersections, especially between many genomes, more clear.

I've always used bbduk.sh to remove host reads. I've used the database linked below, and introduced in the seqanswers post below. One of the benefits to this method is that it separates FASTQ reads into different sets without requiring a bam, so it saves on hard disk space. The method below also uses a masked reference so as not remove true microbial reads that have homology to the human genome. My two comments are: 1. could you benchmark with bbduk as well? I think this could be a real contribution in this space. 2. could you mask the T2T assembly like the bbduk author did? I've really appreciated not accidentally removing plant/fungal etc sequences from my metagenomes.

Seqanswers introduction to method: http://seqanswers.com/forums/archive/index.php/t-42552.html Database link: https://drive.google.com/file/d/0B3llHR93L14wd0pSSnFULUlhcUk/edit?usp=sharing

I think it might be worthwhile to create a simulated set of reads from the human genome and from some microbes, such as E. coli or others with high quality genomes where it can be confirmed that there is no human contamination therein. Especially for ribosomes and other sequences with homology, it could be very difficult to establish a ground truth. The inquiries you did with real data are super important as well, and well designed given that real data is messy, but I think simulated data would be a strong contribution here.

It would be really interesting to see the alignment rates -- e.g. what fraction of each sample assembled, and if this varies by biome. This would give us some sort of idea if there were other viral reads left on the table due to non-assembly

Did you do any sort of contamination screen here to see if any of your hits were off target or had homology to other sequences? Either against BLAST nt/nr or against metagenomes or something?

Did you do validation here? We've recently done something similar and noticed that we have to filter our diamond BLAST-equivalent results to 90% identity, or else we get a ton of off target hits.

clever. Reminds me of NCBI's new clustered nr database

Nice! I have two questions about this. 1. Are there any problems that could arise in training because this training set is so unbalanced? 2. How do your input RdRPs compare to those used in Serratus?

I saw that noarch was specified on conda, but when I tried to install it via conda on an M1 mac, I encountered issues whether running on arm64 or rosetta.

I think this sentence has a missing word

This is interesting. I'm curious if this could also stem from contamination in the databases, as mentioned in the next sentence. Is there a way to systematically evaluate this? (e.g. potential for kit contamination in isolates vs assembly/binning contamination in MAGs)

If I'm reading this correctly, this design wouldn't allow you to estimate overfitting. Have you brainstormed any ways to make a validation set for this model?

Would it be possible to add context into the abstract about whether this is a new database, a combination of existing data bases, or just an existing database? I think that would be useful information to have up front

How does this set of genomes compare to e.g. all bacteria and archaea in GenBank? Was there a reason for excluding GenBank?

The genome-grist tool does something similar but is built around sourmash gather and BWA for mapping: https://github.com/dib-lab/genome-grist

I'm curious what the advantage of mapping with an LCA algorithm is here. We've found LCA methods to lead to higher false positives for genome identification (see here: https://www.biorxiv.org/content/10.1101/2022.01.11.475838v2). It would be really helpful to have this explained a bit more.

This is a really clear explanation!

Would you be willing to provide a reference for this? Maybe https://journals.asm.org/doi/10.1128/mSystems.00020-16 or https://bluegenes.github.io/2022-paper-protein-kmers/manuscript.pdf

What statistics does this refer to?

This matches intuition, but is super useful to have it written out so clearly.

would it be possible to contribute parts of your pipeline back to the nf-core eager pipeline? It think it could give it more visibility and be more reproducibly than your snakemake pipeline given the use of docker containers etc. Your approach seems very nifty and would be great to make it more available! Parts of this might be simplified by work being done on the nf-core taxprofiler workflow which I think is building modules for krakenuniq

I took a look at your repository and it was very easy to navigate and follow! I have a couple of suggestions though. Would you be able to pin versions of software tools recorded in your envs/*yaml folder? This will make your pipeline more reproducible.

It would also be helpful if you could parameterize your output directory. A lot of clusters separate write directories from run directories, so this would help with portability.

Lastly, I'm curious if k-mer trimming would improve your results. khmer's trim-low-abund.py can trim samples with variable coverage using the -V flag if you're interested in exploring that further!

I'm so excited your code is open source! Your website is so beautiful I assumed the code would be hidden. Very excited to explore everything there.

Would it be possible to make it so that test data sets have pre-calculated results? I wanted to demo your web service, but the results took a very long time to generate even though I used the Demo Fileset CRC_health.100

Also is it possible for the query task status page to auto refresh the status? it stated running until I refreshed and saw that it was finished.

Is this the readily-measurable environment? Could there be things that are predictive of the environment that weren't measured but impact evolutionary trajectories?

Are the environmental regimes truly identical, or highly similar? The use of identical here makes me think of highly controlled laboratory settings instead of different locations

When I follow the links detailed in Belser et al, I can only find coral amplicon sequencing and aerobiomes, plus two other targeted sequencing experiments. Would you be willing to provide the exact accessions for this data in this paper? It would make the data much more findable!

In the last paragraph of the introduction, the authors state that samples were taken from 33 sites. Is it 33 or 32?

It would be super helpful if you could report the fraction of each metagenome that aligned to each coral species in the results or as a supplement

if this was based on metatranscriptomes, what percent of the contigs were from non-eukaryotic species? was percent were from non-coral species?

Could you use a synonym here? This acronym makes it difficult to dive directly into the discussion and interpret the most salient points

Many of your figures have both red and green in the palette, would you be willing to change to a color-blind friendly palette to make your figures more accessible?

This data set is truly such a gift to the scientific community. I love how this paragraph gives a large scale overview of all of the data collected.

This zenodo is under restricted access, would you be willing to open access now that the preprint is posted?

Is this because these corals were sampled side-by-side? If that's the case, I think stating that more explicitly in the abstract, intro, and here would be really helpful signposting to the readers (and would negate my initial comment on "identical" in the abstract)

I'm not sure if this is helpful, but especially if you end up with specific genomes that you want to look for, you could try using sourmash branchwater: https://www.biorxiv.org/content/10.1101/2022.11.02.514947v1. If you have a eukaryotic genome you're interested in, you could sketch it (sourmash sketch) and then use the branchwater tool to search most metagenomes in the SRA to see which ones have high containment with the genome your searched. You could then use the SRA metadata tables to filter to wastewater samples and the dig in more to the biogeography of those.

This would be really helpful context to have in the introduction, since it would inform why you chose to structure the methods (short kb contigs) the way you did.

Would you be willing to briefly describe the size of k-mer used for this? I could imagine very different results for k-mer size of 4 (tetranucleotide abundances) vs. 21 or 31 (which are generally genus or species specific)

These are such exciting data sets! I tried searching for them on the SRA and the ENA the accessions could not be found. Would you be willing make your data public now that you have released your preprint?

I have only worked with two single cell data sets so please forgive me if my comment is naive, but given that you suggest that NOMAD+ works on raw single cell reads and the proclivity of raw sequencing data to contain microbial contamination, do you think it would be beneficial to recommend that users screen their single cell data for microbial contamination prior to running NOMAD+? One rapid way to do this would be with the sourmash gather algorithm (https://www.biorxiv.org/content/10.1101/2022.01.11.475838v2). Again, I'm not sure how contamination in single cell compares to that of bulk RNA seq, but this is an interesting manuscript that profiles what is in the left over fractions of transcriptomes that don't map to the human genome: https://genomebiology.biomedcentral.com/articles/10.1186/s13059-018-1403-7

missing punctation

This sentence was a little hard to understand, I found myself reading it a few times to make sure I grasped the meaning...but it's an important point! I wonder if it can be rephrased? maybe removing "from" at the beginning of the sentence would help?

Is there any information that does exist that would allow you to have a better understanding of this? For example, I'm guessing a lot of the microbes you looked at have other isolate genomes sequences in public databases; could you potentially do a pangenome analysis to help understand the trends you identified? (not suggesting you do this for this paper! But I think this is interesting to brainstorm around and might warrant an additional sentence of potential future methods someone could use to answer this type of question)

This is a really helpful statement and something I was really curious about -- contamination is always a super hard thing to separate from HGT. I wonder if it would be easy to filter this analysis to only long read genomes, and to look at the identity of nearby proteins to the candidate HGTs. Barring misassemblies/chimeric contigs, it might be easier to trace the HGT events in longer reads. Probably not necessary here, just brainstorming about alternate signals that could really solidify that these candidates are HGT and not contam!

I found myself a little confused by the box in the second row on the left hand side. Is this an outline of results as well? Or did you determine which 582 plant associated proteins to screen? or is this 582 bacteria?

This is a little hard to interpret as written. Does this suggest that Arabidopsis DET2 was the acceptor of Leifsonia DET2, or that the transfer may have occurred at the LCA between Arabidopsis and rice and barley?

Additional details about the meaning of the labels would be helpful. I was confused especially by "3. HGT between eukaryotes and PA bacteria"

typo?

Very helpful statement! I think this also says something interesting about your methods. For "ancient" HGT events, it doesn't matter that much what plant genome you look at to identify them

typo?

Very cool/clever!

It could be interesting to complement this type of analysis by looking at GC content or tetramernucleotide frequency in this gene across genomes, and how that compares to host genome background. How similar are these genes to the rest of the genes in a given bacteria's genome?

This is fairly stringent phred score trimming. I'm curious why you used such a high threshold. I'd also be curious to see if dropping the phred score to ~10 increased viral genome recall in these samples

Although this would be a separate effort from what you have reported here, I think it could be fascinating to expand this type of inquiry beyond the GEOTRACES data set. Using a tool like sourmash branchwater, you could use the viral genome database you describe above and search for metagenomes on the sequence read archive that contain those genomes. It could be a very nice complement to see what sort of global diversity/biogeography there is for these types of viruses. https://www.biorxiv.org/content/10.1101/2022.11.02.514947v1

I'm curious if you've ever explored alpha diversity estimation metrics that don't require subsampling. The StatDivLab out of UW has developed some really lovely software for doing this type of thing. DivNet in particular might work as a drop in replacement to subsampling, allowing the full amount of data captured in each sequencing run to be used while accounting for unequal sequencing depth: https://github.com/adw96/DivNet

This is still fairly confusing and I'm not sure what it means

Did you do anything to deal with contamination? Contamination is fairly widespread, even in refseq genomes, and might lead to unexpected results.

shouldn't this be a part of the results section?

can you explain what you mean by novel? not found by other studies? only found in one genome?

This might not be low enough to detect orthology if the HGT event is ancient. One recent paper showed similarity drops as low as ~40% https://doi.org/10.1101/2022.08.25.505314

This paragraph was a bit confusing to follow but I think I got the gist of it after a few passes through! I'm curious if you thought about controlling for natural variation in 4mer frequency throughout the genome, as some other methods have found that this helps reduce off target predictions (reviewed in https://doi.org/10.1371/journal.pcbi.1004095). It may not be necessary since you do a second step after the initial screen, but I was just curious if that was something you thought about putting in place, and if so, why you decided against it

Would you be willing to provide a more clear definition of non-redundant here? does this mean there are no paralogs of the gene? or the HGT only occurred in one model org? or only one genome of all of the 824+13 that you investigated?

Are all of these long read genomes? If not, will this be an unreliable estimate?

How did you assess these numbers? in metagenome binning, 1kb isn't large enough to get confident estimates of tetramernucleotide frequency; you often need > 2500 bp.

evalue can change with the size of the database, how did you account for this?

Is 1 base pair a large enough overlap to be biologically important?

Does this mean the 4 you found are different from the ones found in 21 and 24? if so, how do you account for missing the ones found in 21 and 24?

Why did you use euclidean distance?

It would be super helpful here to point out where these protein sequences come from -- NCBI nr, mmseqs sets, etc.

Can you add a few more details on how you identified this? why are tools in mmseqs2 not sufficient here? what innovation is needed to overcome whatever barriers exist?

is each of the 1.8 million genomes to come from a separate species, or will some be different strains of the same species?