Introduction

Information processing in the neocortex and hippocampus relies on fast excitatory synaptic transmission between pyramidal neurons. Many excitatory synapses that release glutamate project onto tiny structures protruding from dendritic shafts termed spines, and the morphological properties of these spines strongly influence synaptic transmission and function, including the activity-dependent synaptic plasticity implicated in learning and memory^{1, 2}. In general, larger spines contain a greater number of glutamate receptors and are functionally stronger³. The morphological parameters defining spine shape, such as length, radius, and volume, are also highly dynamic. For instance, small nascent spines can be enlarged by high-frequency stimuli inducing long-term potentiation (LTP) of synaptic transmission^{4, 5}; alternatively, spine volume can be reduced by low-frequency stimuli that induce long-term depression (LTD)^{5, 6}. Due to these dynamic changes, spines on pyramidal neurons exhibit substantial morphological variation, even within the same dendritic branch. In addition to direct synaptic activation inducing LTP or LTD, spine morphology can be altered by cyclical hormonal changes, neuromodulators, and various activity-independent processes^7–10. In turn, these changes in spine morphology are associated with marked changes in neural function and behavior.

Recent advances in super-resolution microscopy have provided new opportunities for the study of spine structural diversity and dynamics. A recent investigation using stimulated emission depletion (STED) imaging reported an increase in spine neck size after LTP, resulting in stronger electrical coupling between the spine head and dendritic shaft¹¹. A structured illumination microscopy (SIM) imaging study of cultured hippocampal neurons revealed expansion of the spine head convexity after LTP induction¹², and suggested that this structural transformation may increase the adhesion of presynaptic and postsynaptic membranes. SIM imaging has also been applied to record the dynamic properties of spinules¹³, thin protrusions extending from existing spines that may function in secondary synapse formation. Collectively, these reports implicate spine nanostructural dynamics in multiple core synaptic properties and functions, including neurotransmission, mechanical stability, synaptogenesis, and plasticity.

Genetic mutations and variants that disrupt synaptic function are strongly implicated in the pathogenesis of neuropsychiatric disorders such as autism spectrum disorder (ASD) and schizophrenia^14–20. Based on these genomic aberrations, numerous mouse models of ASD and schizophrenia have been established and subsequently shown to harbor abnormalities in both spine morphology and physiology^{19, 21–24}. Moreover, several studies have reported synaptic dysfunction in neurons differentiated from patient-derived pluripotent stem cells ^25–28. Schizophrenia and ASD share several clinical features, suggesting common genetic risk factors and etiology^{29, 30}. However, a diametrical model emphasizes the distinct social-cognitive properties of these two major psychiatric disorders^{31, 32}. Both ASD and schizophrenia also show highly heterogeneous clinical phenotypes, potentially stemming from distinct pathogenic processes within disease categories. Therefore, in addition to clinical features, comparative biological phenotyping of multiple ASD and schizophrenia models may reveal shared as well as disease-specific pathogenic processes.

Genetic studies on ASD and schizophrenia have identified copy-number variants and protein-disrupting single-nucleotide variants that confer higher disease risk, and mouse models harboring these genetic abnormalities frequently exhibit behavioral and neurological features resembling those of the clinical conditions (termed endophenotypes)^{19, 33}. However, there are only a limited number of studies objectively comparing phenotypes related to synaptic structure and function between disease models. An in vivo two-photon imaging study reported accelerated spine turnover in three genetically distinct ASD mouse models ³⁴. Also, multiple ASD mouse models were found to exhibit a similar imbalance in the activity of excitatory and inhibitory neurons³⁵. These studies support the idea that distinct genetic mutations may induce similar synaptic and circuit phenotypes in ASD mouse models. However, differences in synaptic properties between ASD and schizophrenia have not yet been comprehensively examined.

In this study, we developed an objective method for identifying population-level differences in spine nanostructure. This method, which has been applied to multiple mouse models of psychiatric disorders, identified two distinct groups corresponding to ASD and schizophrenia. An increase in a specific spine subpopulation with small volumes characterized the spines of schizophrenia mouse models. In turn, ASD mouse models showed a tendency of decreasing spine subpopulation with small volumes. The schizophrenia-related phenotype was associated with slower volume growth and enhanced turnover. Further gene expression analysis identified the overexpression of Ecrg4, a gene encoding a precursor of hormone-like peptides, as a candidate mediator of this schizophrenia-associated spine phenotype. Population-level analysis of spine nanostructures is a powerful approach for understanding heterogeneous synaptic impairments in psychiatric disorders.

Results

Spine nanostructure imaging and generation of the spine density plot

We have designed a method for objectively comparing spine properties among multiple mouse models of psychiatric disorders at the population level using SIM imaging. Previous results demonstrated that high-resolution SIM images of dendritic spines in cultured hippocampal neurons can provide reliable information about spine size and morphological properties¹². In the current study, the analytical procedures for SIM images were further improved for quantitative comparisons. Three-dimensional SIM images of dendritic spines were captured in DiI-stained cultured neurons derived from control and disease-model mouse hippocampus (Figure 1A). This staining procedure provided a better image quality than fluorescent protein-based labeling because of the higher fluorescence signal from the plasma membrane. SIM images were processed using custom-made scripts to identify and segment the individual dendritic spines (Figure 1B–D). We previously reported a set of scripts designed to isolate dendritic spines automatically¹². However, the segmentation accuracy still necessitated manual inspection. In the newly developed scripts, we introduced a step for readjusting the spine–shaft boundary that enabled fully automated measurement of multiple morphological parameters, including length, surface area, and volume. In addition, the isolated spines were further divided into 160-nm thick longitudinal segments for calculating a total of 64 shape descriptors.

Figure 1

Spine shape parameters measured by conventional methods are difficult to compare directly when images are obtained from multiple independent experiments. This problem is caused by several reasons, but the main factors are variability in image intensity, signal-to-noise ratio, and image aberration introduced by unoptimized imaging conditions. To test the stability of the SIM-based spine imaging, we compared the spine size distribution in four samples prepared at time points separated by >2 months from each other (Supplementary Figure 1). The cumulative distribution curves of the four datasets showed extensive overlap, confirming the reproducibility of the SIM-based spine shape measurements using independent culture samples.

Another technical issue related to the comparison of spine shape is the heterogeneity of spine shape and the absence of defined subclasses. Population-level spine analysis, therefore, requires a new strategy that enables objective comparisons of large numbers of spines with a continuum of morphologic features. We previously demonstrated that principal component analysis (PCA) was effective for comparing high-dimensional spine shape features (descriptors) in a feature space with reduced dimensions¹². Indeed, a PCA-based dimensional reduction applied to a large number of spines based on their descriptors generated a two-dimensional plot of spine shape that closely matched conventional spine subtypes, such as thin, mushroom, and stubby (Figure 1E). To objectively compare the two spine datasets (control and disease model), we first generated plots of relative spine densities (number per area in feature space) mapped to 80 × 80 square blocks (Figure 1F), and then subtracted the corresponding plot pairs defined by genotype (control vs. mutant) or culture conditions. The resulting subtracted density plots contain comprehensive information about the differences in relative spine numbers with similar shape properties.

Comparison and grouping of multiple mouse models based on morphological spine similarity

The subtracted density plots contain rich information about the population-level differences in spine shape between samples. Therefore, the generation and comparison of the plots for both ASD and schizophrenia mouse models will be useful for their unbiased grouping. We selected the following eight representative mouse models harboring either copy-number variation or single gene mutations associated with psychiatric disorders: male hemizygous or female homozygous mutant mice for neuroligin-3 R451C mutation (Nlgn3^{R451C/(y or R451C)})^{36, 37}, heterozygous Syngap1 mutant mice (Syngap1^+/−)³⁸, heterozygous POGZ mutant mice (POGZ^Q1038R/+)³⁹, mice paternal duplication of chromosome 7C (15q11-13^dup/+) corresponding to human 15q11-13 duplication⁴⁰, mice with heterozygous deletion of chromosome 16B2,3 (3q29^del/+) corresponding to human 3q29 deletion⁴¹, mice heterozygous for deletion on chromosome 16qA13 (22q11.2^del/+) corresponding to 1.5 Mb deletion on human 22q11.2⁴², mice heterozygous for Setd1a mutation (Setd1a^+/-)⁴³, and Ca²⁺/calmodulin-dependent protein kinase IIα kinase-dead mutant mice (CaMKIIα^K42R/K42R)⁴⁴. Previous studies on these mouse models confirmed the presence of multiple sensory, memory, and social endophenotypes of psychiatric disorders compared to controls. Of these, Nlgn3^{R451C/(y or R451C)}, POGZ^Q1038R/+, and 15q11-13^dup/+ models mimic genetic variations found in ASD. Human mutations in the SYNGAP1 gene are associated with mental retardation (MR) and autistic behaviors. The 3q29^del/+, 22q11.2^del/+, and Setd1a^+/- models harbor genetic mutations associated with higher schizophrenia risk. The CaMKIIα^K42R/K42R mutation does not correspond to any known human genetic variant, but CaMKII-related signaling pathway disruption has been implicated in the working memory deficits found in schizophrenia patients^{45, 46}. Therefore, these mouse models cover a broad spectrum of human genetic aberrations associated with psychiatric disorders. Direct comparison of spine phenotypes could therefore provide important clues to underlying shared and disease-specific synaptic pathologies.

For these comparisons, we prepared three independent culture samples from each control–mutant pair and obtained 700–1,500 3D-SIM spine images 3 weeks after plating for both mutant and control samples (Nlgn3^{R451C/(y or R451C)}, n = 1,134 and 1,204; Syngap1^+/−, n = 991 and 1,371; POGZ^Q1038R/+, n = 1,341 and 1,271; 15q11-13^dup/+, n = 1,208 and 1,099; 3q29^del/+, n = 1,429 and 1,408; 22q11.2^del/+, n = 1,143 and 914; Setd1a^+/-, n = 1,381 and 1,405; CaMKIIα^K42R/K42R; n = 880 and 763, control and mutant, respectively). Initial comparison of the cumulative frequency plots for spine length, spine surface area, and spine volume alone did not provide sufficient information to infer their structural similarity (Supplementary Figure 2), indicating that simple structural parameters are not effective in identifying disease-specific spine shape characteristics. We next generated subtracted density plots from 8 mouse models. These plots were variable, but there were certain similarities between the specific mutants (Figure 2A). We evaluated pairwise similarities by calculating the 2D cross-correlations of the subtracted density plots. The matrix of these 2D cross-correlation values yielded two major groups with high within-group similarity (Figure 2B), and the presence of these two groups was further confirmed by unbiased clustering analysis (Figure 2C). The first group contained the ASD models Nlgn3^{R451C/(y or R451C)}, Syngap1^+/−, and POGZ^Q1038R/+, while the second group consisted of the schizophrenia models 3q29^del/+, 22q11.2^del/+, Setd1a^+/-, and CaMKIIα^K42R/K42R. This result is striking, as it indicates that the information of spine shape alone can tell whether the culture preparation is derived from ASD or schizophrenia mouse models. It is also of note that there were greater similarities among the schizophrenia mouse models than among the ASD mouse models. This finding suggests a convergence of spine pathophysiology in the schizophrenia mouse models.

Figure 2

In addition to the eight mouse models of psychiatric disorders, we generated subtracted density plots for three additional conditions: immature culture samples (13 days in vitro, DIV), treatment with the AMPA subtype glutamate receptor blocker CNQX, and treatment with the GABA_A receptor blocker bicuculline (Figure 2A). The subtracted density plots of the 13 DIV culture and CNQX treatment were similar to those of the schizophrenia mouse models. This finding suggests that the hippocampal synapses in mouse models of schizophrenia are immature and inactive. In contrast, the subtracted density plot of the bicuculline-treated sample was distinct from that of both the schizophrenia and ASD mouse models.

Distinct spine shape between ASD and schizophrenia mouse models

The subtracted density plots provide information about the structural properties of specific spine subsets, which either increased or decreased in disease mouse models. We first defined the areas within the plots where the control or mutant spines showed higher density. Next, we extracted the spines in the control-enriched area (A) and the mutant-enriched area (B), and compared the volumes and shape profiles of spines in these areas. For the 22q11.2^del/+ schizophrenia mouse model, the control-enriched area (A) contained larger spines than the mutant-enriched area (B) (Figure 3B, D, F, H, and J). In contrast, the mutant-enriched area (B) of the Nlgn3^{R451C/(y or R451C)} ASD mouse model overlapped with the domain of large spines, while the control-enriched area (A) contained smaller spines (Figure 3A, C, E, G, and I). These differences in population-level spine properties between ASD and schizophrenia mouse models were relatively consistent (Supplementary Figures 3–7) in that smaller spines were observed in mutant-enriched areas of all four schizophrenia mouse models (3q29^del/+, 22q11.2^del/+, Setd1a^+/-, and CaMKIIα^K42R/K42R) (Supplementary Figures 6), while large spines were dominant in the mutant-enriched areas of the three ASD mouse models (Nlgn3^{R451C/(y or R451C)}, Syngap1^+/−, and POGZ^Q1038R/+). The only exception was 15q11-13^dup/+, which showed enrichment of small spines in the mutant-dominant area. The fractions of spines present in the control-enriched area (A) and mutant-enriched area (B) were substantial in all eight mouse models (fraction of spines in area A: 0.24 ± 0.058 for controls, 0.20 ± 0.063 for mutants; fraction of spines in area B: 0.26 ± 0.057 for controls, 0.31 ± 0.062 for mutants, mean+SD), indicating that despite broad value spectra for all datasets, approximately half of the total spine population was within either area (A) or (B) in data from both the control mice (0.50 ± 0.036, mean+SD) and mutant mice (0.51 ± 0.031, mean+SD).

Figure 3

Slow initial growth of dendritic spines in schizophrenia model mice

The subtracted density plots indicate the common changes in the spine population of the schizophrenia mouse models. The mutant-enriched areas correspond to the specific spine population with small morphology, suggesting underlying dysfunction in either generation of nascent spines or shrinkage of existing spines. To examine potential abnormalities in spine dynamics, we performed time-lapse imaging of hippocampal neuron cultures derived from two schizophrenia mouse models, 22q11.2^del/+ and Setd1a^+/-, and the Nlgn3^{R451C/(y or R451C)} ASD mouse model as a reference (Figure 4). Time-lapse imaging over 24 h captured numerous de novo spine formation and spine loss events. Compared to control neurons, neurons from the two schizophrenia models exhibited slower growth of nascent spines (Figure 5G–L). Spine turnover rate was upregulated in these schizophrenia mouse models (Figure 5A and B), while the ASD model did not show significant upregulation (Figure 5C). The slow growth of spines in neurons from the schizophrenia mouse models was associated with a longer lifetime of transient spines (Figure 5D and E) compared to the ASD model (Figure 5F). Taken together, these results suggest that the slower growth of nascent spines and instability of existing spines contribute to the development of schizophrenia-related spine phenotypes.

Figure 4

Figure 5

To further examine the relationship between the spine growth rate and lifetime, we performed simulations using the computational model shown in Figure 6A. In this model, new spines appear randomly along the dendritic shaft according to a density-dependent probability function, then increase in volume with speed V1 (phase 1), are eliminated with probability P1 (phase 2), and are converted to stable spines upon reaching a volume threshold (phase 3). These stable spines are subsequently destabilized with low probability P2 (phase 4) and start to shrink with speed V2 (phase 5). When these shrinking spines reach a threshold, they are permanently removed. This model was run over 10 days of neuronal differentiation in culture (corresponding to the period of spine density increase from 8 to 17 DIV) with 625 different combinations of V1, V2, P1, and P2 to identify the combination with the best fit to the experimental data for spine turnover of 15–18 DIV neurons derived from mutants 22q11.2^del/+, Setd1a^+/-, and Nlgn3^{R451C/(y or R451C)} over 24 h (Figure 6B). The best-fitted models all yielded spine lifespan distributions matching observations of the three mouse mutants (Figure 6C). The three characteristics of the spines in 22q11.2^del/+ and Setd1a^+/-were also preserved in the best-fitted model. First, the simulation showed higher turnover rates for schizophrenia-related spines (21.4 + 2.89% and 36.2 + 4.16% for the best-fitted models for the control and the schizophrenia mouse models, Supplementary Figure 8A), as in the case of experimental data (Figure 5A and B). Second, spine density in the schizophrenia mouse models was lower in both the simulation and experimental data (simulation; 0.57 + 0.029 /μm and 0.40 + 0.027 /μm for the control and the best-fitted for schizophrenia, shown in Supplementary Figure 8B, experiment; 0.69 + 0.17 /μm and 0.56 + 0.20 /μm for the control and 22q11.2^del/+, 0.70 + 0.15 /μm and 0.64 + 0.16 /μm for the control and Setd1a^+/-). Third, the size of transient spines was smaller in the schizophrenia mouse models in the simulation (0.33 + 0.033 and 0.20 + 0.024 (a.u.) for the control and the best-fitted for schizophrenia, Supplementary Figure 8C), consistent with the data of spine nanostructure (Figure 3H).

Figure 6

After the selection of the best-fitted parameters for three mouse models, we compared the optimal values of V1, V2, P1, and P2 between mouse models (Figure 6D). The nascent spine volume increase (V1) was slower in the schizophrenia models compared to corresponding controls. In turn, the probability of nascent spine loss (P1) was higher in the ASD model. In all three disease models, the probability of spine destabilization (P2) was moderately higher than that in the corresponding controls. These results suggest that a slow spine volume increase in schizophrenia model mice results in the accumulation of nascent spines with smaller volumes. On the other hand, the spines of ASD model mice grow at a rate similar to that of control spines, but are more likely to be eliminated both in the nascent stage and after maturation.

The simulated trajectory plots of individual spines along the dendritic shaft using the best-fit parameters also reflected the differential patterns of spine turnover in culture preparations derived from the control and schizophrenia model mice (Figure 6E). In culture preparations derived from schizophrenia mouse models, there were fewer spines with a lifetime longer than 2 or 3 days, while there was a greater number of spines formed and eliminated within 24 h. This increase in the number of less stable spines explains the pattern of spine lifetime in Figure 5D and E. In summary, these simulations suggest that reduced spine maturation rate can account for the three key properties of schizophrenia-related spines: the greater number of small volume spines, the enhanced turnover, and the reduction in spine density observed in schizophrenia model mice.

Screening of molecules responsible for altered spine properties in schizophrenia mouse models

Both SIM-based imaging and time-lapse imaging suggest that the initial growth phase of spines is slower in schizophrenia model mice. To identify candidate genes that may contribute to this growth impairment, we performed RNA sequencing of libraries derived from the cultured neurons of eight psychiatric disorder mouse models and corresponding controls, and filtered the results for differentially expressed genes (DEGs). The DEGs specific to more than two schizophrenia-related mouse models are listed in Supplementary Figure 9. In the hippocampi of schizophrenia model mice, the expression levels of Met, Osr1, Arhgap15, and Rnf170-ps were downregulated, while Cip4, Ecrg4, and Npas4 expression levels were upregulated. The majority of these DEGs have been previously implicated in the pathology of schizophrenia (Met)⁴⁷, synaptic plasticity (Npas4)^{48, 49}, or signaling pathways in the nervous system (Arhgap15⁵⁰, Cip4⁵¹, and Ecrg4⁵²).

If these genes do regulate spines in schizophrenia mouse models, manipulation of these genes in wild-type neurons may induce a spine phenotype similar to schizophrenia mouse models. To test this possibility, we downregulated Met and Arhgap15, and upregulated Cip4, Ecrg4, and Npas4 in separate wild-type hippocampal neuron cultures and measured spine dynamics using time-lapse imaging. Ecrg4 was the only DEG that modulated spine dynamics in the predicted direction (Figure 7A and B). Thus, the aberrant expression of Ecrg4 may be related to the spine phenotype of schizophrenia model mice.

Figure 7

Ecrg4 encodes a precursor protein that yields several short hormone-like peptides via multiple cleavage sites⁵³. Ecrg4 has been reported to function in cell proliferation, inflammation, and tumor suppression. However, there is limited information on the neuronal functions of Ecrg4. Therefore, we further analyzed the cellular distribution of Ecrg4 in wild-type hippocampal neurons to assess its potential functions at synapses. The Ecrg4 protein was detected in the membrane fraction of cultured hippocampal neurons (Supplementary Figure 10A). In addition, HA-tagged Ecrg4 showed punctate localization patterns within the neuronal cytoplasm (Figure 7C). A fraction of HA-tagged Ecrg4 puncta localized at axonal varicosities (arrowheads) or dendritic spines (arrows), together with some extracellular clusters (asterisks). HA-tagged Ecrg4 also colocalized with GFP-tagged neuropeptide Y (NPY) (Supplementary Figure 10B and C), suggesting transport in dense-core vesicles. Surface labeling of live cells also revealed clusters of Ecrg4 on axonal boutons and dendritic spines (Supplementary Figure 10D and E), supporting active exocytosis and tethering via its binding partners to the plasma membrane. Collectively, these observations suggest that the intracellular machinery for transporting and releasing Ecrg4 may be actively involved in the regulation of synapse formation and maturation.

Rescue of the synaptic phenotype in schizophrenia mouse models by Ecrg4 suppression

The spatial localization of Ecrg4 protein in wild-type neurons and upregulated gene expression in the hippocampus of schizophrenia model mice strongly suggest a contribution to the observed spine abnormalities. To examine this possibility, we tested whether Ecrg4 suppression rescues the synaptic phenotype of neurons derived from schizophrenia mouse models. We introduced a GFP expression plasmid along with the Ecrg4-shRNA plasmid or the control shRNA plasmid into neurons derived from 22q11.2^del/+ and Setd1a^+/- mice (Figure 8A and Supplementary Figure 10G). While the spine density plots generated from GFP-expressing neurons differed from those generated from DiI-labeled neurons, the subtracted (differential) distribution pattern between Ecrg4 shRNA and control shRNA neurons showed an increase in the larger volume spine fraction (Figure 8B), suggesting that Ecrg4 suppression reversed the spine phenotype observed in naïve neurons from schizophrenia model mice. In contrast, the subtracted density plot between control shRNA-expressing neurons from schizophrenia model mice and wild-type mice revealed a reduction in the fraction of large volume spines (Supplementary Figure 11), consistent with comparisons of DiI-labeled naïve neurons from these genotypes. These data support the idea that the shRNA-mediated downregulation of Ecrg4 can rescue the spine shape impairment observed in schizophrenia mouse models.

Figure 8

Discussion

We developed an objective method for identifying population-level differences in spine nanostructures and applied it to cultured hippocampal neurons derived from multiple mouse models of psychiatric disorders. We succeeded in identifying two distinct groups of mouse models according to the population-level spine properties. Notably, these two objectively identified groups correspond to mouse models with ASD-related gene mutations and schizophrenia-related gene mutations. The higher fraction of small spine population was the common characteristics of schizophrenia mouse models. Time-lapse imaging of spine dynamics revealed that this smaller spine phenotype is associated with a reduced speed of volume increase and enhanced turnover. Screening for genes differentially expressed between mutant and control mice identified Ecrg4 overexpression as a potential molecular mechanism contributing to the schizophrenia spine phenotype. Indeed, suppression of Ecrg4 expression normalized the spine phenotype in neurons derived from schizophrenia mouse models. These findings support population-level analysis of spine nanostructures as a powerful approach to elucidating the heterogeneous synaptic impairments underlying psychiatric disorders.

Spine nanostructures provide information about developmental history, previous activity, and current functional properties, but these fine structural details were not accessible by conventional light microscopy. However, the recent advances of super-resolution microscopy allow for the study of nanoscale spine morphology, including in living preparations, and have revealed new structural features related to important synaptic functions^11–13. Another challenge to spine nanoscale imaging is the structural heterogeneity across large spine populations, which may obscure important population-level differences. Nonetheless, to understand the physiology and pathology of neural networks, it is essential to obtain comprehensive information about a large number of spines belonging to a single neuron and within local neuronal networks. We previously reported a dimensional reduction technique coupled with machine learning for evaluating spine properties at the population level ¹², and have further refined this method for the current investigation to extract population-level nanostructural features of spines in multiple mouse models of psychiatric disorders. This method enabled us to compare spine phenotypes among different mouse models objectively and group disease models according to similarities and differences in the spine population. Notably, the two identified groups corresponded to ASD-related gene mutations and schizophrenia-related gene mutations, supporting the idea that population-level analyses of spine nanostructures provide meaningful information on synaptic abnormalities and dysfunction.

Population-level spine properties were more homogeneous in the schizophrenia models (those with gene mutations implicated in schizophrenia) than in the other 4 models studied, due in part to a shared tendency for smaller spines. In contrast to these findings, several previous studies reported reduced numbers of small spines in the postmortem cortical tissues of schizophrenia patients^{22, 54}. However, this discrepancy may be attributed to developmental period as we examined fetal neuronal spines after three weeks of in vitro differentiation while the aforementioned studies included postmortem samples from patients over 40 years of age. One possible scenario that may help unify these observations is that spine density is reduced in the early postnatal period, and that this reduction evokes homeostatic plasticity, which in turn reduces the fraction of small spines in later life. Another potential explanation is that schizophrenia differentially affects regional spine properties (e.g., hippocampus versus cortex). Previous in vivo two-photon imaging studies reported higher spine turnover in the hippocampus than the cortex^{55, 56}. Therefore, the spine phenotype identified in this study may be specific to hippocampal pyramidal neurons. Finally, while the spine phenotype identified in the human postmortem brain undoubtedly resulted from complex interactions among genetic background, environmental influences, and regulation by nonneuronal cells, the data obtained from pure neuronal cultures are more likely to reflect the direct impact of schizophrenia-related gene mutations on synaptic properties. Therefore, our analytical method is advantageous in that it can reveal the direct effects of genotype on spine phenotype, in this case identifying Ecrg4 overexpression as a potential driver of spine abnormalities in schizophrenia.

In addition to Ecrg4, these analyses identified several other DEGs common to mice harboring schizophrenia-related mutations. However, time-lapse imaging identified Ecrg4 as the most effective spine regulator, and we also found that Ecrg4 suppression reversed the abnormal spine phenotype in neurons derived from schizophrenia model mice. Ecrg4 encodes a precursor protein of hormone-like peptides ⁵⁷. Exogenously expressed Ecrg4 exhibited a punctate pattern in both axonal and dendritic compartments, with a fraction of clusters present on the cell surface. This observation indicates that Ecrg4 may be secreted locally by axons or dendrites, thereby influencing synapses and adjacent neurons. Several signaling pathways act downstream of Ecrg4, including the stress-related NF-κB pathway and PI3K/Akt/mTOR pathway⁵³. Ecrg4-derived peptides may also act through Toll-like receptor 4 (TLR4) or multiple scavenger receptors, such as LOX1⁵⁷. While Ecrg4-related signaling pathways have been studied primarily in the context of inflammation and infection, there are few studies on Ecrg4-related signaling pathways in neurons. However, an Ecrg4-deficient mouse model demonstrated enhancements in neural stem cell proliferation and spatial learning, suggesting that overexpression may disrupt hippocampal function⁵². We suggest that Ecrg4 may be a critical regulator of neural development and synaptogenesis in the hippocampus. Further studies are warranted to identify Ecrg4 receptors and downstream signaling pathways and to investigate the potential dysfunction of these pathways in schizophrenia mouse models.

The nanoscale features of dendritic spines in ASD mouse models were more variable than those in schizophrenia mouse models. This difference may be related to the broader clinical spectrum of ASD, which ranges from mild impairments in social skills to severe intellectual disability. It is possible that the four ASD mouse models examined in this study, Nlgn3^{R451C/(y or R451C)}, Syngap1^+/−, POGZ^Q1038R/+, and 15q11-13^dup/+, represent subgroups with different levels of hippocampal dysfunction. Among the four ASD mouse models, 15q11-13^dup/+ showed population-level spine properties closer to those of the schizophrenia models, and several clinical studies have reported that CNVs in the 15q11-13 chromosomal region also confer an increased risk of schizophrenia^{58, 59}. Analysis of the relationships between rare genetic variants and synapse phenotypes in mouse models may contribute to the eventual categorization of diseases based on neural circuit function rather than heterogenous and overlapping symptoms.

The most important advantage of our analytical method is reproducibility, despite the substantial heterogeneity in genotype and phenotype among these mouse models. Thus, we are confident that the spine population analysis described here will help researchers in different laboratories obtain reliable and comprehensive information on spine properties. Further, this method could facilitate comparisons among a large number of mouse models harboring genetic mutations related to mental disorders. Spine morphology as measured in a reduced culture system can also demonstrate direct effects of gene mutations on neuronal phenotypes in the absence of indirect influences from nonneuronal cells or specific environments. Ultimately, these studies may identify the shared and specific pathogenic mechanisms of psychiatric disorders and lead to improved diagnosis and treatments.

Materials and methods

Neuronal culture and genotyping of tissue samples from genetically modified mice

The following mouse models of psychiatric disorders were examined for spine nanostructural characteristics: male hemizygous or female homozygous neuroligin-3 R451C mutant (Nlgn3^{R451C/(y or R451C)}), heterozygous Syngap1 mutant mice (Syngap1^+/−), heterozygous POGZ mutant mice (POGZ^Q1038R/+), mice with duplication of chromosome 7C (15q11-13^dup/+; corresponding to human 15q11-13 duplication), mice heterozygous for deletion on chromosome 16B2,3 (3q29^del/+; corresponding to human 3q29 deletion), mice heterozygous for deletion on chromosome 16qA13 (22q11.2^del/+; corresponding to 1.5 Mb deletion on human 22q11.2), mice heterozygous for Setd1a mutantion (Setd1a^+/-), and Ca²⁺/calmodulin-dependent protein kinase IIα kinase-dead mutant mice (CaMKIIα^K42R/K42R). All animal experiments were reviewed and approved by Institutional Animal Care and Use Committee of the University of Tokyo (permission number A2023M031-02).

For the production of embryos harboring heterozygous mutations, male heterozygous mutants were bred with wild-type female C57BL/6J mice (Japan SLC). For the production of Nlgn3^{R451C/(y or R451C)} embryos, Nlgn3^R451C/y males and Nlgn3^R451C/+ females were crossed. Embryos of CaMKIIα^K42R/K42R mice were obtained by crossing heterozygous CaMKIIα^K42R/+ males and females. Both male and female embryos were used for the primary culture. Control cultures were prepared from embryos without disease-related mutations on the background of C57BL/6J.

Genomic DNA was purified from E16 embryos prior to brain dissection using the QuickGene DNA tissue kit (WAKO), and genotypes were determined by PCR using either KapaTaq (Kapa Biosystems) or KOD FX Neo (TOYOBO) following the standard protocols provided by the manufacturers. The primer sequences used for genotyping were as follows:

Nlgn3_R^{451C/(y or R451C)} genotyping:

Sense primer: 5′- GGTCAGAGCTGTCATTGTCAC-3′

Antisense primer: 5’- TGTACCAGGAATGGGAAGCAG-3’

Syngap1^+/- genotyping:

Sense primer for WT: 5’-GTCAGTGGGACATGGAAGTAG-3’

Sense primer for mutant: 5′-CTTCCTCGTGCTTTACGGTATC-3′

Antisense primer (common): 5’-CTGATCAGCCTGTCAGCAATG-3’

POGZ^Q1038R/+genotyping:

Sense primer for WT: 5’- TCTGTGAAGAAGCTTCGGGTAGTAC-3’

Antisense primer for WT: 5′-GTCTCCTCATTTACAGGGAGCT-3′

Sense primer for mutant: 5’- GCAGCGGCTCCCCGTTAAC-3’

Antisense primer for mutant: 5′- AGCGCACAGCCCACTCATAG-3′

15q11-13^dup/+ genotyping:

Sense primer: 5’- AGAGGAGGGCCTTACTAATTACTTA-3’

Antisense primer: 5′-ATATGTACTTTTGCATATAGTATAC -3′ 3q29^del/+ genotyping:

Sense primer for WT: 5’- TTGGCACCACTCGCCCAAGTTATATCCACC-3’

Sense primer for mutant: 5’- CAAGGGGGAGGATTGGGAAGACAATAGCAG-3’

Antisense primer (common): 5’- GGTCATGCAAATTCTAGCAGTGAGTCATGAC-3’

22q11.2^del/+ genotyping:

Sense primer for WT: 5’- GAGAAAGTGGGTGGGAAGGC -3’

Antisense primer for WT: 5’- GTCCCTCGCCACAGTCATAA -3’

Sense primer for mutant: 5’- CTAAGGAATGGTTCCGGCCA -3’

Antisense primer for mutant: 5′- TTTCACGGAGGCGGTATTCA-3′

Setd1a^+/- genotyping:

Sense primer: 5′-CTCGCCGCCATTTCTCTACATC-3′

Antisense primer: 5′-GTTCTGGAGGTTCTGGAGGTG-3′

CaMKIIα^K42R/K42R genotyping:

Sense primer: 5’-GGTCTTGAAGACTGTCTGGTGTGAGA-3’

Antisense primer: 5′-CACAGGCCAGTTTAGGTCTTGCTAGG-3′

After genotyping, primary hippocampal neuron cultures were prepared from mutant and corresponding control lines as described previously^{12, 60}. In brief, E16 hippocampi were dissected, minced, treated with trypsin (Gibco) and DNase (SIGMA), and dissociated mechanically into a cell suspension. Cells were resuspended in minimum essential medium supplemented with B18 supplement, L-glutamine (Gibco), and 5% fetal calf serum (Equitech-Bio), and plated onto glass-bottom dishes (MatTek, #1.5) precoated with poly-L-lysine (SIGMA). Glial proliferation was prevented by adding 5 μM Ara-C (SIGMA) 2 days after plating. For the comparison of spine nanostructure between wild-type and mutant neurons, culture samples were processed for imaging at 18-22 DIV. For pharmacological treatment, 10 μM CNQX or 20 μM bicuculline was applied to the culture medium at 13 DIV and samples were processed for imaging at 19 DIV.

Transfection and fluorescence labeling of primary hippocampal neurons

Primary neurons were transfected with the indicated vectors using the calcium phosphate precipitation method at 8–9 DIV according to a previously described procedure⁶¹. All fluorescent proteins were expressed under the control of the β-actin promoter.

The expression of the shRNA constructs was induced using the pSilencer vector system. Briefly, the 19-nt target sequence and the reverse complement sequence of the target were ligated into the shRNA expression vector separated by a 9-nt spacer loop (TCTCTTGAA) (pSilencer 2.0-U6; Invitrogen). Stealth RNAi Negative Control Duplexes (Invitrogen) were used as control shRNAs. The targeting sequences for mouse Met, Arhgap15, and Ecrg4 are as follows:

Met RNAi sequence

5’-GCAGTGAATTAGTTCGCTA-3

Arhgap15 RNAi sequence:

5’-AAGACAGATGTGAACATAC-3

Ecrg4 RNAi sequence:

5’- GAGGCTAAATTTGAAGAT-3’

The suppression of target protein expression by these transfected shRNAs was tested by immunoblotting lysates of COS-7 cells co-transfected with a target gene construct fused to a sequence encoding amino acids 98–106 of human influenza hemagglutinin (HA), all under the control of the β-actin promoter. Overexpression of Cip4, Ecrg4, and Npas4, genes upregulated in cultures derived from mouse models of schizophrenia, was induced by transfection of plasmids containing the corresponding HA-fused cDNA under the control of the β-actin promoter. Briefly, mouse Trip10 transcript variant 4 (Cip4; NM_001242391) and mouse Ecrg4 (NM_024283) were cloned, and an HA tag was inserted into the C-terminus. Mouse Npas4 (NM_153553) was cloned according to a previous report⁶² and fused with a Flag-HA tag at its C-terminus. Human neuropeptide Y (NPY) tagged with green fluorescent protein (GFP, a kind gift from Dr. J. Takaska, NIH) was expressed under the control of the β-actin promoter. For cell surface labeling of individual neurons in culture, cells were fixed and stained with 1,1’-dioctadecyl-3,3,3’,3’-tetramethylindocarbocyanine perchlorate (DiI; Molecular Probes) as described previously¹².

Immunocytochemistry

Cultured hippocampal neurons expressing HA-tagged constructs were fixed with 2% paraformaldehyde in PBS, permeabilized with 0.2% Triton X-100 in PBS, blocked with 5% normal goat serum, and reacted with rabbit anti-HA antibody (1:500; MBL) and mouse anti-β-galactosidase antibody (1:500, SIGMA). After washing with PBS, samples were reacted with Alexa Fluor 546-labeled goat anti-rabbit IgG antibody (1:500; Invitrogen) and Alexa Fluor 647-labeled goat anti-mouse IgG antibody (1:500; Invitrogen). For cell surface labeling with GFP, neurons expressing Ecrg4-HA-SEP (superecliptic pHluorin) were reacted with anti-GFP VHH/nanobody conjugated to Alexa Fluor 647 (1:4,000; Chromotek) in culture medium at 37°C.

Structured illumination microscopy (SIM) imaging

For 3D-SIM imaging, we used an inverted microscope (ECLIPSE Ti-E, NIKON) equipped with an N-SIM-E system, 405, 488, 561, and 640 nm diode lasers (LU-NV, NIKON), and an oil-immersion TIRF objective lens (SR Apo TIRF 100 ×, N.A. 1.49, NIKON). All SIM images were acquired under strictly controlled conditions, including a stage temperature of 28–29°C, to minimize the position and aberration fluctuations. Spherical aberrations induced by refractive index mismatch were reduced by manually adjusting the objective lens correction collar. Images were acquired using an EMCCD camera (iXon3 DU-897E, Andor Technology) with 512 × 512 pixels (each pixel, 16 μm square) operated in the read-out mode at 1 MHz with 16-bit analog-to-digital conversion of the EM gain. All image-processing steps were performed in three dimensions. The acquired image datasets were computationally reconstructed using the reconstruction stack algorithm V2.10 of NIS-Elements AR (NIKON). An image stack consisted of 63 axial (z) planes of 512 × 512 pixels in the x- and y-dimensions, with pixel size of 32 × 32 nm. All 3D images were reconstructed in 120-nm z-steps spanning 7.56 μm in the z-axis. The voxel size fulfilled the Nyquist criterion requirements. Only dendritic segments spatially isolated from other fluorescent objects were selected for analysis because objects brighter than the targeted dendritic segments made image thresholding unreliable. Quality checks of the acquired SIM images were performed as described previously using the SIMcheck plugin for ImageJ⁶³.

Spine isolation and shape measurement

Spines were isolated by automated image thresholding using a previously described method¹². Briefly, the reconstructed SIM image stacks were first processed by multilevel thresholding using the MATLAB built-in function [multithresh()] to produce binary images, which were further processed by geodesic active contours using the MATLAB built-in function [activecontour()]. The binary image stacks generated by thresholding the SIM images were then processed for the automated detection of spines. First, the dendritic shafts were fitted with elliptic cylinders and voxel clusters outside the best-fit elliptic cylinders were identified as spine candidates. Next, the spine candidates were classified by their volumes and shapes. After selection, the junction between the dendritic shaft and spine was determined. Finally, the nanostructural parameters of isolated spines were measured. These steps are basically the same as those in our previous algorithm, but several points were improved, as described in the following sections (The original code for SIM image processing is available at https://shigeookabe.github.io/download-page-SIM/ with a password of “simspineimage”).

(1) Dendritic shaft fitting

We first estimated the shaft volume using the custom-made function “SIM_Fitdendrite” as previously described. In the new algorithm, the candidate spines were subsequently evaluated to determine whether their bases were large and spread on the shaft surface. If the spine base was large, the basal voxels were removed using the function “SIM_Removevoxel.”

(2) Selection of the spines

Using the custom-made function “SIM_Position,” spines with volumes too small, too elongated, or located at the image edge were removed as previously described. Small spine candidates with volumes less than twice the volume of a rectangular cuboid with edge lengths equal to the theoretical resolutions of SIM images (<0.01 μm³) were rejected.

(3) Adjustment of the spine–shaft junction

We determined the location of spine–shaft junction using the custom-made function “SIM_Junction.” This step was identical to our previous algorithm. Next, spines with more than one connection site to the shaft were detected using the function “SIM_Neck” and eliminated from further analysis.

(4) Measurement of nanostructural parameters

First, the angle of the spines was adjusted to be perpendicular to the long axis of the dendritic shaft and within the xy-plane using the function “SIM_Angle.” Next, nanostructural parameters were measured using the function “SIM_Shape.” In our previous spine shape measurement study, a polygon mesh-based calculation was performed. In the current study, we took a different approach using voxel-based measurement. The calculation of structural parameters using voxel data was more straightforward and required fewer calculation steps than the method based on polygon mesh data. Spines were divided into 160-nm segments along the long axis using the function “SIM_Shape.” These spine segments were used to calculate the following 64 nanostructural parameters: spine length, surface area, base surface area, total volume, volume of each spine segment (20 segments), convex hull volume of each spine segment (20 segments), and convex hull ratio of each spine segment (20 segments).

Analysis of the spine population data

Differences in spine phenotype among model mice were analyzed at the population level by PCA. We selected five parameters (spine length, base surface, total volume, and volumes of the fifth and fifteenth spine segments) for the PCA-based dimensional reduction. Two parameters that reflect the principal structural features (length and volume) were first selected. Second, three other parameters that were mutually independent and also independent from the first two parameters were chosen (pairwise correlation coefficients < 0.3). The first three principal components (PC1–PC3) covered approximately 90% of the variance in the data.

The 2D presentation of spine population data using PC1 and PC2 reflected the overall pattern of the spine shape distribution, with thin, mushroom, and stubby spines separated in the feature space. Next, relative spine densities (number per area in the feature space) were mapped onto 80 × 80 square blocks spanning the standard deviation of ±4. Subsequently, two density plots, one derived from culture samples of disease model mice and the other from corresponding control mice, were subtracted to extract the differences in spine numbers with similar shape properties. The similarity between the two subtracted density plots was evaluated by calculating the 2D cross-correlation C using the following equation

where X and Y correspond to the two subtracted density matrices. The 11 different subtracted density plots shown in Figure 2 were clustered using the MATLAB function “linkage,” and the average distance between all pairs of elements in the two samples was calculated using the cluster algorithm.

Areas in the feature space enriched with mutant or control spines were then extracted by searching the 80 × 80 square blocks and identifying those where the difference in relative density was larger than the mean plus 3 × SD. Spine lengths, volumes, and axial profiles were then analyzed within areas enriched with either mutant or control spines.

Time-lapse imaging of dendritic spines and image analysis

Live cell 3D imaging was performed using an A1 confocal laser scanning microscope (NIKON) equipped with 405, 488, 561, and 640 nm diode lasers (LU-NV, NIKON). Cells at 15–18 DIV were maintained at 37°C under a continuous flow of 5% CO₂ to maintain medium pH using a heater stage system (INUG2H-TIZSH, Tokai Hit) with a lid designed for a glass-bottom dish to minimize evaporation. Horizontal and vertical drifts were controlled by a motorized stage and a perfect focus system, respectively. Images were acquired with an oil immersion objective lens (60×, NA) 1.4, NIKON) at a confocal aperture of 1 AU. The imaged volume was 30.7 × 30.7 × 7.56 μm in x, y, and z directions, respectively, and image stacks consisted of 17 planes separated by 500 nm. Time-lapse images were recorded either at short intervals (15 min intervals for 24 h) or long intervals (two time points separated by 24 h). Cultured neurons imaged at 24-h interval were then fixed and immunostained with anti-HA to confirm the expression levels of the transfected genes.

The image stacks obtained by time-lapse imaging were analyzed using custom-made MATLAB scripts. First, 2D images were generated by maximal intensity projection of z-stack images. These 2D images were processed to transform the outlines of dendritic shafts to be aligned with straight lines. Spines were subsequently identified as objects extending beyond the straight edges of the dendritic shafts. Identified spines were segmented, and total fluorescence intensities were measured as a proxy for spine volume. We classified the spines as stable (present throughout the imaging period), newly formed (appearing within 24 h imaging period), and eliminated (lost within 24 h imaging period). For newly formed and eliminated spines, we generated temporal profiles of growth and shrinkage using the total fluorescence intensity values.

RNA sequencing and gene expression analysis

Total RNA was isolated from cultured hippocampal neurons at 19 DIV using the RNeasy Mini Kit (QIAGEN) with the RNase-Free DNase Set (QIAGEN) according to the manufacturer’s protocols. The quality of the extracted RNA was initially assessed using a NanoDrop spectrophotometer, and the RNA integrity number (RIN) score was calculated for each sample using a Bioanalyzer (Agilent Technologies). The yields were quantified using a Qubit Fluorometer (Thermo Fisher Scientific). Sequencing libraries were prepared using the NEBNext Ultra Directional RNA Library Prep Kit for directional libraries (New England BioLabs) and the KAPA HTP Library Preparation Kit (KAPA Biosystems) according to the manufacturers’ protocols. Paired-end RNA-seq libraries were sequenced using the Illumina HiSeq and MGI DNBSEQ-G400 platforms. Three independent culture preparations were used for RNA sequencing of each mouse model and control library. Following quality control, the reads were aligned with the mouse reference genome GRCm38/mm10 (GenBank assembly ID: GCA_000001635.2). Expression levels of genes were analyzed based on the transcripts per million and differential expression analysis was conducted using the DESeq2 package in R⁶⁴. Differentially expressed genes (DEGs) were selected according to an adjusted P < 0.05 and a (log2)-fold change greater than 0.5.

Simulation of spine development

Spine growth was simulated using a dendrite model consisting of 1,000 segments, each 0.25 μm in length (250 μm in total). The spine distribution was updated every 20 min over 10 days, covering the initial phase of spine development from 8 DIV to 17 DIV. The initial probability of new spine generation was set at 7.5 × 10⁻³ per segment per hour and decreased thereafter in proportion to the accrued spine density increase (by the final 24 h, the probability of spine generation [mean ± SD] was 3.2 × 10⁻³ ± 0.23 × 10⁻³/segment/h). This probability function yielded an average of 1.2 ± 0.22 new spines per 10 μm of dendrite every 24 h, in good agreement with the spine dynamics measured by time-lapse imaging from 15–18 DIV (1.3 ± 0.63 new spines per 10 μm of dendrite every 24 h). The model also set the increase in spine volume (V1) from 6.7 × 10⁻³ to 3.3 × 10⁻² in 20 min. Spines in the early growth phase were treated as unstable, with loss probability P1 ranging from 9.0 × 10⁻³ to 0.1 every 20 min. Spines were modeled as stabilized when the individual volume reached a relative value of 1.0. These stabilized spines then entered a shrinking phase with a probability P2 ranging from 2.9 × 10⁻³ to 4.0 × 10⁻³ in 20 min. The shrinking rate V2 was set from 6.7 × 10⁻³ to 3.3 × 10⁻² per 20 min. In simulations, the four parameters V1, V2, P1, and P2 were varied within these indicated ranges to find the combination producing the best fit with experimental data on the difference in spine lifetime distribution profile and turnover rate at 24 h (the actual measurement interval). The simulation with the best-fit parameters yielded a spine density comparable to measured data (0.57 ± 0.029 vs. 0.70 ± 0.17 spines /μm of dendrite).

Western blotting

Total protein extracts from cultured hippocampal neurons were separated into cytosolic and membrane fractions using Trident Membrane Protein Extraction Kit (GeneTex) according to the manufacturer’s protocol and separated on 10% sodium dodecyl sulfate (SDS)-polyacrylamide gels or 15%–20% Tricine-SDS-polyacrylamide gels. The separated proteins were then transferred onto nitrocellulose membranes (Millipore) or PVDF membranes (Bio-Rad) using a wet transfer system. Membranes were blocked with 5% bovine serum albumin in Tris-buffered saline plus Tween-20 (TBS-T, 20 mM Tris-HCl, 137 mM NaCl, 0.2% Tween-20) for 1 h at room temperature, then incubated overnight at 4°C with rabbit anti-Ecrg4 (1:1,000, SIGMA), rabbit anti-GluA1 (1:1,000, Nittobo Medical), and mouse anti-α-tubulin (1:1,000, SIGMA) as indicated. The membranes were washed and incubated with HRP-conjugated secondary antibody against mouse or rabbit IgG (1:10,000, Amersham) for 2 h. Signals from the target protein bands were detected using a SuperSignal™ West Atto detection kit (Pierce) and imaged for quantitation (densitometry) using a ChemiDoc imaging system (Bio-Rad).

Statistics

The statistical tests used for each experiment and the sample numbers (of neurons, dendrites, or spines) are indicated in the corresponding figure legends. Statistical significance was determined using the Wilcoxon rank sum test for two independent samples (Figures 3, 5, and 8), two-way ANOVA for genotype and time (Figure 5), or one-way ANOVA with post hoc Dunnett’s tests between groups (Figure 7). Statistical significance was set at P < 0.05.

Resource availability

The original code for SIM image processing is available at https://shigeookabe.github.io/download-page-SIM/ with a password of “simspineimage”.

Supplementary Figures

Supplementary Figure 1.

Supplementary Figure 2.

Supplementary Figure 3.

Supplementary Figure 4.

Supplementary Figure 5.

Supplementary Figure 6.

Supplementary Figure 7.

Supplementary Figure 8.

Supplementary Figure 9.

Supplementary Figure 10.

Supplementary Figure 11.

Acknowledgements

We thank Toru Takumi for the 15q11-13^dup/+ mice, Katsuhiko Tabuchi for Nlgn3^{R451C/(y or R451C)} mice, Yoko Yamagata for CaMKIIα^K42R/K42R mice, and Justin Taraska for the NPY tagged with GFP. This work was supported by Grants-in-Aid for Scientific Research (20H00481, 20H05894, 20H05895 to S.O.), the Japan Agency for Medical Research and Development (JP19gm1310003 to S.O., T. N., and Y. G., JP22jm0210097 to S.O., JP19dm0207071 to A. A.), the Naito Foundation, and the Uehara Memorial Foundation.

Additional information

Author contributions

Y.K. and S.O. conceived of the project and designed the methodology. Y.K. and Q.L. performed the experiments. R. S., A.A. and T.N. provided the animal resources and provided advice on analysis. Y.G. conducted the gene expression analysis. S.O. and Y.K. wrote the manuscript, with contributions from A.A., T.N., and Y.G.

Declarations of generative AI and AI-assisted technologies in the writing process

The authors utilized Grammarly (https://app.grammarly.com/) for editing the text. The authors take full responsibility for the contents of the manuscript.