Does it make biological sense to interpret peptide length as a potential source of systematic difficulty in ACE-inhibition prediction (different presentation / models (regression) ) or should I avoid making that interpretation and frame this purely as a data/prediction phenomenon??? and what if i cross veify those peptides on different publicly available datasets ?

Hi all,

I'm looking to use GROMACS to perform a molecular dynamics simulation of a poly beta amino ester (PBAE) polymer assembling with its mRNA cargo in order to obtain values like radius of gyration and hydrophobic and hydrophilic surface area.

I know that GROMACS requires a .pdb file to start with, but how could I get that for a PBAE polymer? I'm also not sure how to go about doing this in GROMACS (I'm pretty new to molecular dynamics simulations) so any guidance would be helpful!

i can help myself with a blog or something

I want to investigate whether a transcriptionally distinct intermediate state exists between two previously characterized cell populations.

The relevant scRNA-seq dataset only has raw SRA reads available..no processed matrix or Seurat object. Reconstructing it would require ~150 GB of FASTQs and a custom reference, which could take significant time and compute.

My main concern is that there is only one library per condition, so the reconstructed data wouldn't provide biological replication for strong condition-level statistical comparisons.

However, cell-level data may be necessary for my specific question, since the available cluster summaries aren't sufficient.

Would reconstructing the dataset be worthwhile, or is there a sensible way to scope the analysis down while still answering the biological question?

Hi everyone. I am new to bioinformatics in general. I am a biochem currently doing a bioengineering phd (still in pre-candidature). I had some snippets of bioinformatics during my undergrad but nothing beyond BLAST and docking. Never had formal programming formation, just side projects and AI-guided R coding for small data analysis and graphs.

For what I want to do for my thesis I really need to learn omics analysis properly, specially transcriptomics. During self-learning, I stumbled upon this amazing resource (https://www.sc-best-practices.org/) on single cell transcriptomics, so I have been following it as my starting point and learning cool stuff, thank you to the authors of it!

Anyways, since I've already had some experience with R, I decided to try and learn python bioinformatics as an excuse to learn python too. In the interoperability section of the book I mentioned the authors state

"A common question from new analysts is which ecosystem to focus on (referring to Bioconductor, Seurat or the Scverse). While it makes sense to start with one, and a successful analysis can be performed in any ecosystem, competent analysts should be familiar with all three and comfortable moving between them. This allows analysts to always use the best-performing tools, regardless of their implementation. Analysts who are not comfortable switching ecosystems often default to familiar packages, even when better alternatives exist elsewhere"

Which makes sense and sounds logical good advice. But then, doing exercises on public GEO datasets on bulk RNA-seq as practice, still with the mindset of sticking to python as an excuse to learn it, i stumbled upon an article (Colange et al. 2025 here) of a project that migrates a lot of tools of R to the scverse. In there, authors rationale is that python is the new default language everyone learns and they create the library InMoose to migrate or directly replace, for example, DESeq2. Furthermore, besides direct drop-in replacement tools, the authors frame python as the future choice (at least, as part of the rationale).

So, as a guy who is just starting, I wanted to ask people with experience in bioinformatics (you all) either developers or tool-users:

1) Do you marry an ecosystem like scverse or Bioconductor and just work in there for comfort? Or do you switch frequently depending on the needs?

2) For people who doesn't come from an informatics background, how long did it take for you to learn your niche and what were your best resources/helpers?

3) Do you think python will ever replace R in data analysis?

4) What is your opinion on AI-guided learning? (as for me, I use gemini to solve questions or create graphics presets but sometimes by seeing other people's codes I realize that it mashes up some concepts or methods from various pipelines into a coherent-resulting graph that I am not always sure if they make sense)

5) Do you create your own pipelines/portfolio to analyze data? Or you just tweak existing pipelines?

I appreciate any answer to any of those questions, thanks for your time in at least reading

Hi- I'm really stuck on the best way to proceed with snpEff for annotating my VCFs. Prior to running snpEff, I created an "ancestral" fasta based on several outgroups and polarized my VCFs to that fasta. But now, when building the snpEff database, I'm not sure if I should be using that ancestral fasta or my original reference genome. I tried both ways and unfortunately am getting different counts of variant types. So, it clearly makes a difference, but I'm not sure which way is correct. Thanks!

Hi all,

I’ve been working with a Xenium 500-plex spatial dataset for the past year and, although I’ve been able to make good use of it overall, there are still a few aspects of the analysis that I haven’t fully resolved.

Annotation:
My cell-type annotation is based mainly on manual annotation of clusters from the integrated dataset, as I have not been able to fully annotate the tissue using Tangram because the available gold-standard scRNA-seq references do not capture all of the cell populations present in my samples.

The clusters are generated using graph-based Leiden clustering on Harmony dimensions derived from SCTransform-normalised data. One issue I consistently encounter is that, within most annotated clusters, a proportion of cells falls away from the main cluster centroid on the UMAP. I suspect this may partly reflect spatial transcript bleeding or contamination from neighbouring cells, although spatial bleeding correction itself require good enough references and/or alters counts making use of normalisation like SCTransform not feasible (and my samples have variable quality). My tissue also does not have a well-defined organisation of cell types into discrete anatomical structures, so approaches such as SpaNorm do not seem particularly well suited to the dataset.

DEG analysis of co-localised cells:
A related issue arises when performing differential expression analysis on spatially co-localised cells. For example, when comparing a cell type of interest X depending on whether it is located close to cell type Y, I find it difficult to determine whether increased expression of certain genes in X reflects a genuine transcriptional response of X to its proximity to Y, or simply transcripts originating from Y being incorrectly assigned to neighbouring X cells.

I’d be very interested to hear how others working with Xenium data have approached either of these issues.

I have gone deep into the world of agentic AI and have been applying it to my genomics work. Academic research lab scRNAseq, spatial omics, epigenomics etc I have used all kinds of agent harnesses, skills, connectors and extensions. Everyone seems convinced it is the future of our field and if we don’t get ahead of it we will be left behind. But I feel like my progress is much slower now, much messier I spend a lot of time doing directory archaeology. I understand less of what I produce. I thoroughly don’t enjoy my work anymore. Should I risk getting left behind and give up with AI and agentic bioinformatics, for the sake of my sanity?

I'm relatively new to bioinformatics but interested in doing an independent research project. I have some basic knowledge and am willing to learn the necessary tools and concepts along the way.

How realistic is it to take on a meaningful project independently without much prior experience? What would be a good way to choose a project that’s challenging but still manageable, and what should I learn before starting?

Would appreciate advice from anyone who has done something similar!

I’m building a wheat regulatory state prediction model. There is this CNN model called DeepWheat and I would like to know if it’s possible to turn it into a transformer through transfer learning using my own multiomics tracks and a graph pangenome, as DeepWheat maps their tracks to AK58 and that is not good enough for what I’m building. How would you approach this?

Hi everyone,

I am a PhD student still learning the bioinformatics side of single-cell analysis. I am working with nine independent patient-derived 10x Genomics Multiome samples. I'm reaching you guys because I need some sanity check because I've been going back and forth on this for over two months now.

Since I have 10x Multiome multisamples from the same tissue that respond differently to various treatments, this is the workflow I have used:

After per-sample QC (both RNA and ATAC), after that I merged the object and then:

1. Recalled the peaks separately for each patient with MACS3 (defaults parameters).
2. Combine the nine per-patient peak sets using GenomicRanges::reduce().
3. Keep peaks between 20 and 10,000 bp.
4. Remove non-standard chromosomes and ENCODE-blacklisted regions.

This is what I got first time:

A senior bioinformatician I've been consulting keeps telling me this is wrong, and that it should look more like the second image (green one below, this is not the actual consensus peak) instead — unimodal, right-skewed, with peaks less than the -d MACS3 value, mode around 200–300 bp, with the gap between the median and the mean being roughly 2–3x, and a final consensus peak count somewhere between 100k and 250k due to my number of samples.

I've spent two months trying to find where exactly this recommendation comes from and I haven't been able to track it down clearly in any paper or tutorials or any official page that deals with sc/snATAC-seq. I've asked him repeatedly and keep getting the same answer, so I'm fairly confident I'm not misunderstanding what's being asked of me...

After many weeks of trying, I was able to have an approach that reduces the peak count gradually while moving the width distribution toward the requested range following this exact workflow:
MACS3 command (per sample):
macs3 callpeak -t <fragments.tsv.gz> -f FRAG -g hs -n <sample> --outdir <dir> \

-q 0.05 --min-length 10 --max-count 1 --barcodes <barcodes.txt> \

--call-summits --buffer-size 10000

My unfiltered nine-patient union contains approximately 500k peaks, with:

• median width: ~244 bp;
• mean width: ~430 bp;
• modal bin: below 100 bp

Filtering step using −log10(qvalue) from narrowpeaks.

Here an actual example of my sweep:

There is precedent for filtering or ranking MACS peaks before constructing a union, although the implementations are not identical to mine:

• Lyu et al. ranked peaks according to their −log10(q-value), retained the top 120,000 peaks for each cluster in each sample, and then merged them into the final union peak sets using GenomicRanges::reduce().

So per-sample score/q-value filtering isn't unheard of in the sc/snATAC-seq literature. But I know that in my case I'm using it to reverse-engineer a target shape rather than for a principled reason (FDR, reproducibility, etc.), and that's the part that doesn't sit right with me methodologically.

So, two things I'd really appreciate feedback on:

1. Are these specific numbers (mode 200–300 bp, median/mean ratio ~2–3x) an actual documented convention that I'm just failing to find, or does this sound more like a personal heuristic being presented as a hard rule (he previously analyse the data, so maybe he is expecting that I follow his way(?) ?

2. Is filtering per-sample peaks by q-value before the union/merge step a methodologically defensible way to build a consensus peak set, or is there a more standard approach I should be using instead?

Any papers, pipelines or examples showing the expected consensus peak-width distribution, rather than the fragment-length distribution, would be extremely helpful.

Thanks for your time and consideration!

Hello,

I am interested in comparing protein structures using foldseek. I can do this using query/subject protein structures predicted from sequence (which I have seen other people have done with e.g. ProstT5) and/or I can use 'real' structures from something like PDB.

My concern with real structures is that proteins are dynamic, so I suspect one protein might have multiple different structures in PDB. I could then have a situation where a protein has a conformation (in real life) that is a very good match for my query protein, but that conformation is not stored for that protein (in PDB).

My thinking, then, is that I might be better to use just predicted structures and ignore real structures, because the predicted structures are derived from a uniform source (sequence) and so will have only one predicted confirmation. If I predict everything using the same structure prediction model, I would then expect structurally-similar proteins to have similar predicted structures (even if the predicted structures themselves aren't perfectly accurate). Whereas for PDB, I might have essentially the same protein relative to my query, but different conformations stored for that protein, and those conformations might not be similar to my query at all.

Does that make sense?

Hello all, I want to analyse the interaction between a homotrimer and heterodimer. There are no computational studies on it, but people have biologically validated the interaction by doing a mutagenesis study. My job is to do a computational analysis of the same. I docked them using HADDOCK, giving appropriate restraints from experiments, but the distance between the proposed active sites is around 25+ Angstorms. How do I go further? Is this even a feasible study? I tried docking the heterodimer with 50-kmer, it didn't give me the usual interaction distance (3 Angstorms or so).Help me, I need the expertise of a structural biologist, How do I go further? 😭😭😭

Is there any video or site that explains the E value ?

How did you guys start molecular docking? Any articles or videos that helps you figure everything out? I didn't have any teacher or resources in my city sadly so I have to rely on free resources 😭 Thank you in advance!!

I am trying to assign taxonomy to an 18S amplicon dataset using dada2::assignTaxonomy() with the full PR2 reference database in Rstudio.

Relevant details:

• 147 samples
• ~18,900 ASVs
• PR2 v5.1.1 DADA2-formatted SSU reference
• minBoot = 80
• multithread = 1
• 16 GB physical RAM
• Windows 11, 64-bit R

With the full dataset, RAM reaches 100% and Windows starts heavily paging to disk.

To check whether the number of query ASVs was the main issue, I repeated the taxonomy assignment using only 500 ASVs. RAM still reached ~100% while DADA2 was processing the reference FASTA, followed by very high disk activity from paging.

Has anyone successfully run assignTaxonomy() with the full recent PR2 database?

In particular, I would like to know:

• Is 32 GB RAM usually sufficient?
• Would 64 GB RAM be a safer requirement?
• What peak RAM usage have others observed with recent PR2 versions?

I am mainly trying to determine the appropriate memory requirement before moving the taxonomy step to a larger workstation or cloud instance.

I am thankfull for any kind of advise or nudge in the right direction!

There is a protein that I have to work on for my research. I picked it from PDB and it had many missing residues. So I fixed it with modeller. I have further refined it's Ramachandra plot, errat, molprobity value and some other parameters that are never for validation. But I am having a problem with its 3d structure. It's verify 3d value is not really good so I want to refine it. I have tried some tools but they didn't work. Can you guys please help me with that. If there are any good tools that I can use for it and they should be easy to use as well

Hi everyone, for those who regularly perform scRNA sequencing analysis, how do you all find cell type specific marker genes?

I have been trying to find cell type specific markers from all my clusters in my scRNA data with sc.tl.rank_genes_groups (a scanpy function), which performs a Wilcoxon test to find the genes that are dominantly expressed in one cluster compared to the other clusters, however, I’ve noticed that some of the genes provided in the output are actually highly expressed in more than one cluster, especially among small and similar cluster populations. I think this may due to the averaging effect of the method, which makes it great for cell type populations that are drastically different compared to those that aren’t. Hence, I was wondering if there’s a more accurate method of finding genes that are specific to a single cell type.

Thank you in advance for the advice!

I am concerned with some long-read sequencing. I noticed that certain samples had very large numbers of SVs (particularly inversions). I usually filter reads for average base quality and read length and do QC with cramino. I noticed that median read length is very low (500bp) and that the % mapped reads is low (~60%) but the median percent identity is fine (98%) and coverage is decent (34x). Looking at the read alignments on IGV there are lots of reads with supplementary alignments which might explain the low % read maps. I can see the coverage looks pretty "patchy". Im not sure if I should reject these samples or that these are genuine detections of structural rearrangements ?

I'm new and idk how I make sure that my result is accurate. Do I need to do an extra step? I use autodock vina and I believe there's a box you can set the size. Does the size and position of the box affect the docking result? Is there anything I should keep in mind to make sure that it's right? Thank you in advance!

Hi guys, im new to molecular dynamic world. I wanted to ask few questions about MD. The story begins where i did 1000 ns nsteps for my aptamer-protein. Fast forward to backmapping part, where i wanted to create tracjectory and aa files, I couldnt create some of the files for certain frames.

1. So im not sure it is okay for me to continue the analysis part such as PLIP, TTClust and etc if i dont have these files?

2. And 1 more thing, what if i have all the analysis but the result is not quite satisfy, do i need to repeat the same things with different sequence or can i modify or optimize my aptamer-protein interaction in order to get a better analysis result?

Please help me, thanks!

Hi all!

I'm working with a set of metagenome-assembled genomes (MAGs) and trying to annotate metabolic pathways, I specifically focused on methanogenesis and methane oxidation. I ran two different tools on the same dataset:

1. KEGGdecoder – I used Prodigal for gene prediction, then ran the predicted proteins through KoFamScan, retained only high-confidence hits (marked with *), and used that as input for KEGGdecoder.
2. METABOLIC – This tool uses its own integrated pipeline (Prodigal + HMM-based searches) to predict genes and assign functions directly.

The results are quite different. KEGGdecoder identified complete methanogenesis pathways (e.g., M00357) in several bins, while METABOLIC only found partial steps (e.g., M00357+01) in the same bins, and often didn't detect the full module at all. In some cases, METABOLIC didn't find any methanogenesis-related genes in bins where KEGGdecoder reported a complete pathway.

I also noticed that METABOLIC found genes for methane oxidation (mmoB, pmoABC) in some bins, but KEGGdecoder didn't report those pathways — likely because I filtered out lower-confidence hits before running KEGGdecoder.

My questions:

1. Why would KEGGdecoder find a complete pathway when METABOLIC only finds partial steps or nothing at all in the same bin? Is this due to differences in HMM profiles, cutoffs, or how "module completeness" is defined?
2. METABOLIC seems to detect partial pathways and "potential" functions (using a relaxed approach), while KEGGdecoder appears more conservative (requiring all steps). Is it fair to say that KEGGdecoder is better for confirming pathways at the genome/MAG level, while METABOLIC is more useful for community-level trends?
3. Would it be reasonable to use both tools synergistically — METABOLIC for broad community functional profiling and KEGGdecoder for high-confidence pathway confirmation in individual MAGs — or is there a better recommended approach?

I'm not trying to say one tool is better than the other — I understand they use different philosophies. I just want to understand how to interpret and combine these results correctly for my downstream analysis.

Thanks in advance for any insights!

I'm new to bioinformatics and I have NO idea what affinity, distance from rmsd l.b, and best mode rmsd u.b means. I also have no idea what mode is. Can anyone explain? Thank you in advance