Multi-view data offer advantages over single-view data for characterizing individuals, which is crucial in precision medicine toward personalized prevention, diagnosis, or treatment follow-up. Here, we develop a network-guided multi-view clustering framework named netMUG to identify actionable subgroups of individuals. This pipeline first adopts sparse multiple canonical correlation analysis to select multi-view features possibly informed by extraneous data, which are then used to construct individual-specific networks (ISNs). Finally, the individual subtypes are automatically derived by hierarchical clustering on these network representations. We applied netMUG to a dataset containing genomic data and facial images to obtain BMI-informed multi-view strata and showed how it could be used for a refined obesity characterization. Benchmark analysis of netMUG on synthetic data with known strata of individuals indicated its superior performance compared with both baseline and benchmark methods for multi-view clustering. In addition, the real-data analysis revealed subgroups strongly linked to BMI and genetic and facial determinants of these classes. NetMUG provides a powerful strategy, exploiting individual-specific networks to identify meaningful and actionable strata. Moreover, the implementation is easy to generalize to accommodate heterogeneous data sources or highlight data structures.

Pancreatic intraepithelial neoplasia (PanIN) is a precursor to pancreatic cancer and represents a critical opportunity for cancer interception. However, the number, size, shape, and connectivity of PanINs in human pancreatic tissue samples are largely unknown. In this study, we quantitatively assessed human PanINs using CODA, a novel machine-learning pipeline for 3D image analysis that generates quantifiable models of large pieces of human pancreas with single-cell resolution. Using a cohort of 38 large slabs of grossly normal human pancreas from surgical resection specimens, we identified striking multifocality of PanINs, with a mean burden of 13 spatially separate PanINs per cm 3 of sampled tissue. Extrapolating this burden to the entire pancreas suggested a median of approximately 1000 PanINs in an entire pancreas. In order to better understand the clonal relationships within and between PanINs, we developed a pipeline for CODA-guided multi-region genomic analysis of PanINs, including targeted and whole exome sequencing. Multi-region assessment of 37 PanINs from eight additional human pancreatic tissue slabs revealed that almost all PanINs contained hotspot mutations in the oncogene KRAS , but no gene other than KRAS was altered in more than 20% of the analyzed PanINs. PanINs contained a mean of 13 somatic mutations per region when analyzed by whole exome sequencing. The majority of analyzed PanINs originated from independent clonal events, with distinct somatic mutation profiles between PanINs in the same tissue slab. A subset of the analyzed PanINs contained multiple KRAS mutations, suggesting a polyclonal origin even in PanINs that are contiguous by rigorous 3D assessment. This study leverages a novel 3D genomic mapping approach to describe, for the first time, the spatial and genetic multifocality of human PanINs, providing important insights into the initiation and progression of pancreatic neoplasia.

Oncogene amplification is a major driver of cancer pathogenesis. Breakage fusion bridge (BFB) cycles, like extrachromosomal DNA (ecDNA), can lead to high copy numbers of oncogenes, but their impact on intratumoral heterogeneity, treatment response, and patient survival are not well understood due to difficulty in detecting them by DNA sequencing. We describe a novel algorithm that detects and reconstructs BFB amplifications using optical genome maps (OGMs), called OM2BFB. OM2BFB showed high precision (>93%) and recall (92%) in detecting BFB amplifications in cancer cell lines, PDX models and primary tumors. OM-based comparisons demonstrated that short-read BFB detection using our AmpliconSuite (AS) toolkit also achieved high precision, albeit with reduced sensitivity. We detected 371 BFB events using whole genome sequences from 2,557 primary tumors and cancer lines. BFB amplifications were preferentially found in cervical, head and neck, lung, and esophageal cancers, but rarely in brain cancers. BFB amplified genes show lower variance of gene expression, with fewer options for regulatory rewiring relative to ecDNA amplified genes. BFB positive (BFB (+)) tumors showed reduced heterogeneity of amplicon structures, and delayed onset of resistance, relative to ecDNA(+) tumors. EcDNA and BFB amplifications represent contrasting mechanisms to increase the copy numbers of oncogene with markedly different characteristics that suggest different routes for intervention.

The discovery of subtypes is pivotal for disease diagnosis and targeted therapy, considering the diverse responses of different cells or patients to specific treatments. Exploring the heterogeneity within disease or cell states provides insights into disease progression mechanisms and cell differentiation. The advent of high-throughput technologies has enabled the generation and analysis of various molecular data types, such as single-cell RNA-seq, proteomic, and imaging datasets, at large scales. While presenting opportunities for subtype discovery, these datasets pose challenges in finding relevant signatures due to their high dimensionality. Feature selection, a crucial step in the analysis pipeline, involves choosing signatures that reduce the feature size for more efficient downstream computational analysis. Numerous existing methods focus on selecting signatures that differentiate known diseases or cell states, yet they often fall short in identifying features that preserve heterogeneity and reveal subtypes. To identify features that can capture the diversity within each class while also maintaining the discrimination of known disease states, we employed deep metric learning-based feature embedding to conduct a detailed exploration of the statistical properties of features essential in preserving heterogeneity. Our analysis revealed that features with a significant difference in interquartile range (IQR) between classes possess crucial subtype information. Guided by this insight, we developed a robust statistical method, termed PHet (Preserving Heterogeneity) that performs iterative subsampling differential analysis of IQR and Fisher's method between classes, identifying a minimal set of heterogeneity-preserving discriminative features to optimize subtype clustering quality. Validation using public single-cell RNA-seq and microarray datasets showcased PHet's effectiveness in preserving sample heterogeneity while maintaining discrimination of known disease/cell states, surpassing the performance of previous outlier-based methods. Furthermore, analysis of a single-cell RNA-seq dataset from mouse tracheal epithelial cells revealed, through PHet-based features, the presence of two distinct basal cell subtypes undergoing differentiation toward a luminal secretory phenotype. Notably, one of these subtypes exhibited high expression of BPIFA1. Interestingly, previous studies have linked BPIFA1 secretion to the emergence of secretory cells during mucociliary differentiation of airway epithelial cells. PHet successfully pinpointed the basal cell subtype associated with this phenomenon, a distinction that pre-annotated markers and dispersion-based features failed to make due to their admixed feature expression profiles. These findings underscore the potential of our method to deepen our understanding of the mechanisms underlying diseases and cell differentiation and contribute significantly to personalized medicine.

The venetoclax BCL2 inhibitor in combination with hypomethylating agents represents a cornerstone of induction therapy for older AML patients, unfit for intensive chemotherapy. Like other targeted therapies, venetoclax-based therapies suffer from innate and acquired resistance. While several mechanisms of resistance have been identified, the heterogeneity of resistance mechanism across patient populations is poorly understood. Here we utilized integrative analysis of transcriptomic and ex-vivo drug response data in AML patients to identify four transcriptionally distinct VEN resistant clusters (VR_C1-4), with distinct phenotypic, genetic and drug response patterns. VR_C1 was characterized by enrichment for differentiated monocytic- and cDC-like blasts, transcriptional activation of PI3K-AKT-mTOR signaling axis, and energy metabolism pathways. They showed sensitivity to mTOR and CDK inhibition. VR_C2 was enriched for NRAS mutations and associated with distinctive transcriptional suppression of HOX expression. VR_C3 was characterized by enrichment for TP53 mutations and higher infiltration by cytotoxic T cells. This cluster showed transcriptional expression of erythroid markers, suggesting tumor cells mimicking erythroid differentiation, activation of JAK-STAT signaling, and sensitivity to JAK inhibition, which in a subset of cases synergized with venetoclax. VR_C4 shared transcriptional similarities with venetoclax-sensitive patients, with modest over-expression of interferon signaling. They were also characterized by high rates of DNMT3A mutations. Finally, we projected venetoclax-resistance states onto single cells profiled from a patient who relapsed under venetoclax therapy capturing multiple resistance states in the tumor and shifts in their abundance under venetoclax selection, suggesting that single tumors may consist of cells mimicking multiple VR_Cs contributing to intra-tumor heterogeneity. Taken together, our results provide a strategy to evaluate inter- and intra-tumor heterogeneity of venetoclax resistance mechanisms and provide insights into approaches to navigate further management of patients who failed therapy with BCL2 inhibitors.

Mutant selection windows (MSWs), the range of drug concentrations that select for drug-resistant mutants, have long been used as a model for predicting drug resistance and designing optimal dosing strategies in infectious disease. The canonical MSW model offers comparisons between two subtypes at a time: drug-sensitive and drug-resistant. In contrast, the fitness landscape model with N alleles, which maps genotype to fitness, allows comparisons between N genotypes simultaneously, but does not encode continuous drug response data. In clinical settings, there may be a wide range of drug concentrations selecting for a variety of genotypes. Therefore, there is a need for a more robust model of the pathogen response to therapy to predict resistance and design new therapeutic approaches. Fitness seascapes, which model genotype-by-environment interactions, permit multiple MSW comparisons simultaneously by encoding genotype-specific dose-response data. By comparing dose-response curves, one can visualize the range of drug concentrations where one genotype is selected over another. In this work, we show how N-allele fitness seascapes allow for N*2N-1 unique MSW comparisons. In spatial drug diffusion models, we demonstrate how fitness seascapes reveal spatially heterogeneous MSWs, extending the MSW model to more accurately reflect the selection fo drug resistant genotypes. Furthermore, we find that the spatial structure of MSWs shapes the evolution of drug resistance in an agent-based model. Our work highlights the importance and utility of considering dose-dependent fitness seascapes in evolutionary medicine.

Cancerous tumors may contain billions of cells including distinct malignant clones and nonmalignant cell types. Clarifying the evolutionary histories, prevalence, and defining molecular features of these cells is essential for improving clinical outcomes, since intratumoral heterogeneity provides fuel for acquired resistance to targeted therapies. Here we present a statistically motivated strategy for deconstructing intratumoral heterogeneity through multiomic and multiscale analysis of serial tumor sections (MOMA). By combining deep sampling of IDH-mutant astrocytomas with integrative analysis of single-nucleotide variants, copy-number variants, and gene expression, we reconstruct and validate the phylogenies, spatial distributions, and transcriptional profiles of distinct malignant clones, which are not observed in normal human brain samples. Importantly, by genotyping nuclei analyzed by single-nucleus RNA-seq for truncal mutations identified from bulk tumor sections, we show that commonly used algorithms for inferring malignancy from single-cell transcriptomes may be inaccurate. Furthermore, we demonstrate how correlating gene expression with tumor purity in bulk samples provides the same information as differential expression analysis of malignant versus nonmalignant cells and use this approach to identify a core set of genes that is consistently expressed by astrocytoma truncal clones, including AKR1C3, whose expression is associated with poor outcomes in several types of cancer. In summary, MOMA provides a robust and flexible strategy for precisely deconstructing intratumoral heterogeneity in clinical specimens and clarifying the molecular profiles of distinct cellular populations in any kind of solid tumor.

Striatal projection neurons (SPNs) are traditionally segregated into two subpopulations expressing dopamine (DA) D1-like or D2-like receptors. However, this dichotomy is challenged by recent evidence. Functional and expression studies raise important questions: do SPNs co-express different DA receptors, and do these differences reflect unique striatal spatial distributions and expression profiles? Using RNAscope in mouse striatum, we report heterogenous SPN subpopulations distributed across dorsal-ventral and rostral-caudal axes. SPN subpopulations co-express multiple DA receptors, including D1 and D2 (D1/2R) and D1 and D3. Our integrative approach using single-nuclei multi-omics analyses provides a simple consensus to describe SPNs across diverse datasets, connecting it to complementary spatial mapping. Combining RNAscope and multi-omics shows D1/2R SPNs further separate into distinct subtypes according to spatial organization and conserved marker genes. Each SPN cell type contributes uniquely to genetic risk for neuropsychiatric diseases. Our results bridge anatomy and transcriptomics to offer new understandings of striatal neuron heterogeneity.

The antiviral DNA cytosine deaminase APOBEC3B has been implicated as a source of mutation in many different cancers. Despite over 10 years of work, a causal relationship has yet to be established between APOBEC3B and any stage of carcinogenesis. Here we report a murine model that expresses tumor-like levels of human APOBEC3B after Cre-mediated recombination. Animals appear to develop normally with full-body expression of APOBEC3B. However, adult males manifest infertility and older animals of both sexes show accelerated rates of tumorigenesis (mostly lymphomas or hepatocellular carcinomas). Interestingly, primary tumors also show overt heterogeneity, and a subset spreads to secondary sites. Both primary and metastatic tumors exhibit increased frequencies of C-to-T mutations in TC dinucleotide motifs consistent with the established biochemical activity of APOBEC3B. Elevated levels of structural variation and insertion-deletion mutations also accumulate in these tumors. Together, these studies provide the first cause-and-effect demonstration that human APOBEC3B is an oncoprotein capable of causing a wide range of genetic changes and driving tumor formation in vivo .

Spatial cellular heterogeneity contributes to differential drug responses in a tumor lesion and potential therapeutic resistance. Recent emerging spatial technologies such as CosMx SMI, MERSCOPE, and Xenium delineate the spatial gene expression patterns at the single cell resolution. This provides unprecedented opportunities to identify spatially localized cellular resistance and to optimize the treatment for individual patients. In this work, we present a graph-based domain adaptation model, SpaRx, to reveal the heterogeneity of spatial cellular response to drugs. SpaRx transfers the knowledge from pharmacogenomics profiles to single-cell spatial transcriptomics data, through hybrid learning with dynamic adversarial adaption. Comprehensive benchmarking demonstrates the superior and robust performance of SpaRx at different dropout rates, noise levels, and transcriptomics coverage. Further application of SpaRx to the state-of-art single-cell spatial transcriptomics data reveals that tumor cells in different locations of a tumor lesion present heterogenous sensitivity or resistance to drugs. Moreover, resistant tumor cells interact with themselves or the surrounding constituents to form an ecosystem for drug resistance. Collectively, SpaRx characterizes the spatial therapeutic variability, unveils the molecular mechanisms underpinning drug resistance, and identifies personalized drug targets and effective drug combinations.

Human influenza virus evolves to escape neutralization by polyclonal antibodies. However, we have a limited understanding of how the antigenic effects of viral mutations vary across the human population, and how this heterogeneity affects virus evolution. Here we use deep mutational scanning to map how mutations to the hemagglutinin (HA) proteins of the A/Hong Kong/45/2019 (H3N2) and A/Perth/16/2009 (H3N2) strains affect neutralization by serum from individuals of a variety of ages. The effects of HA mutations on serum neutralization differ across age groups in ways that can be partially rationalized in terms of exposure histories. Mutations that fixed in influenza variants after 2020 cause the greatest escape from sera from younger individuals. Overall, these results demonstrate that influenza faces distinct antigenic selection regimes from different age groups, and suggest approaches to understand how this heterogeneous selection shapes viral evolution.

Tumor-infiltrating macrophages support critical steps in tumor progression, and their accumulation in the tumor microenvironment (TME) is associated with adverse outcomes and therapeutic resistance across human cancers. In the TME, macrophages adopt diverse phenotypic alterations, giving rise to heterogeneous immune activation states and induction of cell cycle. While the transcriptional profiles of these activation states are well-annotated across human cancers, the underlying signals that regulate macrophage heterogeneity and accumulation remain incompletely understood. Here, we leveraged a novel ex vivo organotypic TME (oTME) model of breast cancer, in vivo murine models, and human samples to map the determinants of functional heterogeneity of TME macrophages. We identified a subset of F4/80highSca-1+ self-renewing macrophages maintained by type-I interferon (IFN) signaling and requiring physical contact with cancer-associated fibroblasts. We discovered that the contact-dependent self-renewal of TME macrophages is mediated via Notch4, and its inhibition abrogated tumor growth of breast and ovarian carcinomas in vivo, as well as lung dissemination in a PDX model of triple-negative breast cancer (TNBC). Through spatial multi-omic profiling of protein markers and transcriptomes, we found that the localization of macrophages further dictates functionally distinct but reversible phenotypes, regardless of their ontogeny. Whereas immune-stimulatory macrophages (CD11C+CD86+) populated the tumor epithelial nests, the stroma-associated macrophages (SAMs) were proliferative, immunosuppressive (Sca-1+CD206+PD-L1+), resistant to CSF-1R depletion, and associated with worse patient outcomes. Notably, following cessation of CSF-1R depletion, macrophages rebounded primarily to the SAM phenotype, which was associated with accelerated growth of mammary tumors. Our work reveals the spatial determinants of macrophage heterogeneity in breast cancer and highlights the disruption of macrophage self-renewal as a potential new therapeutic strategy.

Androgen receptor- (AR-) indifference is a mechanism of resistance to hormonal therapy in prostate cancer (PC). Here we demonstrate that the HOX/CUT transcription factor ONECUT2 (OC2) activates resistance through multiple drivers associated with adenocarcinoma, stem-like and neuroendocrine (NE) variants. Direct OC2 targets include the glucocorticoid receptor and the NE splicing factor SRRM4, among others. OC2 regulates gene expression by promoter binding, enhancement of chromatin accessibility, and formation of novel super-enhancers. OC2 also activates glucuronidation genes that irreversibly disable androgen, thereby evoking phenotypic heterogeneity indirectly by hormone depletion. Pharmacologic inhibition of OC2 suppresses lineage plasticity reprogramming induced by the AR signaling inhibitor enzalutamide. These results demonstrate that OC2 activation promotes a range of drug resistance mechanisms associated with treatment-emergent lineage variation in PC. Our findings support enhanced efforts to therapeutically target this protein as a means of suppressing treatment-resistant disease.

In the liver, mitochondria are exposed to different concentrations of nutrients due to their spatial positioning across the periportal (PP) and pericentral (PC) axis. How these mitochondria sense and integrate these signals to respond and maintain homeostasis is not known. Here, we combined intravital microscopy, spatial proteomics, and functional assessment to investigate mitochondrial heterogeneity in the context of liver zonation. We found that PP and PC mitochondria are morphologically and functionally distinct; beta-oxidation was elevated in PP regions, while lipid synthesis was predominant in the PC mitochondria. In addition, comparative phosphoproteomics revealed spatially distinct patterns of mitochondrial composition and potential regulation via phosphorylation. Acute pharmacological modulation of nutrient sensing through AMPK and mTOR shifted mitochondrial phenotypes in the PP and PC regions, linking nutrient gradients across the lobule and mitochondrial heterogeneity. This study highlights the role of protein phosphorylation in mitochondrial structure, function, and overall homeostasis in hepatic metabolic zonation. These findings have important implications for liver physiology and disease.

Immunoglobulin A (IgA) is secreted into breast milk and is critical to both protecting against enteric pathogens and shaping the infant intestinal microbiota. The efficacy of breast milk-derived maternal IgA (BrmIgA) is dependent upon its specificity, however heterogeneity in BrmIgA binding ability to the infant microbiota is not known. Using a flow cytometric array, we analyzed the reactivity of BrmIgA against bacteria common to the infant microbiota and discovered substantial heterogeneity between all donors, independent of preterm or term delivery. We also observed intra-donor variability in the BrmIgA response to closely related bacterial isolates. Conversely, longitudinal analysis showed that the anti-bacterial BrmIgA reactivity was relatively stable through time, even between sequential infants, indicating that mammary gland IgA responses are durable. Together, our study demonstrates that the anti-bacterial BrmIgA reactivity displays inter-individual heterogeneity but intra-individual stability. These findings have important implications for how breast milk shapes the development of the infant microbiota and protects against Necrotizing Enterocolitis.

Neuroblastoma is a pediatric cancer arising from the developing sympathoadrenal lineage with complex inter- and intra-tumoral heterogeneity. To chart this complexity, we generated a comprehensive cell atlas of 55 neuroblastoma patient tumors, collected from two pediatric cancer institutions, spanning a range of clinical, genetic, and histologic features. Our atlas combines single-cell/nucleus RNA-seq (sc/scRNA-seq), bulk RNA-seq, whole exome sequencing, DNA methylation profiling, spatial transcriptomics, and two spatial proteomic methods. Sc/snRNA-seq revealed three malignant cell states with features of sympathoadrenal lineage development. All of the neuroblastomas had malignant cells that resembled sympathoblasts and the more differentiated adrenergic cells. A subset of tumors had malignant cells in a mesenchymal cell state with molecular features of Schwann cell precursors. DNA methylation profiles defined four groupings of patients, which differ in the degree of malignant cell heterogeneity and clinical outcomes. Using spatial proteomics, we found that neuroblastomas are spatially compartmentalized, with malignant tumor cells sequestered away from immune cells. Finally, we identify spatially restricted signaling patterns in immune cells from spatial transcriptomics. To facilitate the visualization and analysis of our atlas as a resource for further research in neuroblastoma, single cell, and spatial-omics, all data are shared through the Human Tumor Atlas Network Data Commons at www.humantumoratlas.org.

Intratumoral heterogeneity (ITH)-defined as genetic and cellular diversity within a tumor-is linked to failure of immunotherapy and an inferior anti-tumor immune response. The underlying mechanism of this association is unknown. To address this question, we modeled heterogeneous tumors comprised of a pro-inflammatory ("hot") and an immunosuppressive ("cold") tumor population, labeled with YFP and RFP tags respectively to enable precise spatial tracking. The resulting mixed-population tumors exhibited distinct regions comprised of YFP+ (hot) cells, RFP+ (cold) cells, or a mixture. We found that tumor regions occupied by hot tumor cells (YFP+) harbored more total T cells and a higher frequency of Th1 cells and IFNγ+ CD8 T cells compared to regions occupied by cold tumor cells (RFP+), whereas immunosuppressive macrophages showed the opposite spatial pattern. We identified the chemokine CX3CL1, produced at higher levels by our cold tumors, as a mediator of intratumoral macrophage accumulation, particularly immunosuppressive CD206Hi macrophages. Furthermore, we examined the response of heterogeneous tumors to a therapeutic combination of PD-1 blockade and CD40 agonist on a region-by-region basis. While the combination successfully increases Th1 abundance in "cold" tumor regions, it fails to bring overall T cell activity to the same level as seen in "hot" regions. The presence of the "cold" cells thus ultimately leads to a failure of the therapy to induce tumor rejection. Collectively, our results demonstrate that the organization of heterogeneous tumor cells has a profound impact on directing the spatial organization and function of tumor-infiltrating immune cells as well as on responses to immunotherapy.

Complete characterization of the genetic effects on gene expression is needed to elucidate tissue biology and the etiology of complex traits. Here, we analyzed 2,344 subcutaneous adipose tissue samples and identified 34K conditionally distinct expression quantitative trait locus (eQTL) signals in 18K genes. Over half of eQTL genes exhibited at least two eQTL signals. Compared to primary signals, non-primary signals had lower effect sizes, lower minor allele frequencies, and less promoter enrichment; they corresponded to genes with higher heritability and higher tolerance for loss of function. Colocalization of eQTL with conditionally distinct genome-wide association study signals for 28 cardiometabolic traits identified 3,605 eQTL signals for 1,861 genes. Inclusion of non-primary eQTL signals increased colocalized signals by 46%. Among 30 genes with ≥2 pairs of colocalized signals, 21 showed a mediating gene dosage effect on the trait. Thus, expanded eQTL identification reveals more mechanisms underlying complex traits and improves understanding of the complexity of gene expression regulation.

COVID-19 has become a global pandemic. Immune dysregulation has been implicated, but immune responses remain poorly understood. We analyzed 71 COVID-19 patients compared to recovered and healthy subjects using high dimensional cytometry. Integrated analysis of ~200 immune and >30 clinical features revealed activation of T cell and B cell subsets, but only in some patients. A subgroup of patients had T cell activation characteristic of acute viral infection and plasmablast responses could reach >30% of circulating B cells. However, another subgroup had lymphocyte activation comparable to uninfected subjects. Stable versus dynamic immunological signatures were identified and linked to trajectories of disease severity change. These analyses identified three "immunotypes" associated with poor clinical trajectories versus improving health. These immunotypes may have implications for therapeutics and vaccines.

Type VIIb secretion systems (T7SSb) in Gram-positive bacteria facilitate physiology, interbacterial competition, and/or virulence via EssC ATPase-driven secretion of small ɑ-helical proteins and toxins. Recently, we characterized T7SSb in group B Streptococcus (GBS), a leading cause of infection in newborns and immunocompromised adults. GBS T7SS comprises four subtypes based on variation in the C-terminus of EssC and the repertoire of downstream effectors; however, the intra-species diversity of GBS T7SS and impact on GBS-host interactions remains unknown. Bioinformatic analysis indicates that GBS T7SS loci encode subtype-specific putative effectors, which have low inter-species and inter-subtype homology but contain similar domains/motifs and therefore may serve similar functions. We further identify orphaned GBS WXG100 proteins. Functionally, we show that GBS T7SS subtype I and III strains secrete EsxA in vitro and that in subtype I strain CJB111, esxA1 appears to be differentially transcribed from the T7SS operon. Further, we observe subtype-specific effects of GBS T7SS on host colonization, as subtype I but not subtype III T7SS promotes GBS vaginal persistence. Finally, we observe that T7SS subtypes I and II are the predominant subtypes in clinical GBS isolates. This study highlights the potential impact of T7SS heterogeneity on host-GBS interactions.