A Multivariate Age-Course Framework for Senescence-Related Gene Expression Profiling in a Mouse Model of Osteogenesis Imperfecta

Abstract

Introduction
Osteogenesis imperfecta is primarily characterised by skeletal fragility, but collagen alterations may also affect extraskeletal tissues, making multitissue ageing-related molecular profiling biologically relevant. Ageing is accompanied by molecular changes that do not occur in isolation, but emerge as coordinated shifts across biological pathways, tissues and time. In experimental studies, however, the common practice of analysing molecular data marker by marker may overlook coordinated expression programmes and tissue-specific developmental trajectories. This issue is particularly relevant in small age-course studies, where complex experimental designs require methods able to preserve the structure of the data while providing interpretable summaries of multigene variation.

Objectives
This study aimed to evaluate whether senescence-related gene expression showed coordinated age-course trajectories in BrtlIV mice carrying the heterozygous Col1a1^G349C/+ mutation, a model of Osteogenesis Imperfecta, compared with wild-type Col1a1^+/+ mice. Specifically, we evaluated the joint expression profile of P16, P21, P53 and LMNB1 across age, genotype and tissue, with the goal of identifying tissue-specific molecular trajectories rather than isolated gene-level differences.

Methods
Expression levels of P16, P21, P53 and LMNB1 were analysed in lung, bone and heart samples collected at 6, 12 and 18 months from wild-type and heterozygous mice. Gene-expression values were log10(x+1)-transformed and z-standardised to obtain comparable marker profiles. The four genes were then analysed as a multivariate molecular signature. Differences in coordinated expression patterns were assessed using permutational multivariate analysis of variance based on distance matrices, with permutations constrained within experimental batch identifiers to account for potential batch-related variability. A global model including age, genotype, tissue and their interactions was complemented by tissue-specific models. Homogeneity of multivariate dispersion and sensitivity analyses using an alternative distance metric were used to assess robustness. Ordination-based visualisation, heatmaps and age-course plots were used to support interpretation. Gene-specific mixed-effects models were fitted as complementary analyses to identify the markers contributing to the multivariate signal.

Results
The global multivariate model showed that senescence-marker expression differed across experimental conditions, with age, genotype and tissue acting jointly rather than independently. The model explained approximately 38% of the total multivariate variation, indicating a substantial reorganisation of the combined gene-expression profile. Tissue-specific analyses revealed that these trajectories were not uniform across organs. Lung and heart showed clearer genotype-dependent age-course patterns, whereas bone displayed more heterogeneous expression dynamics. Higher-order interaction terms supported the presence of tissue-specific age-related differences between wild-type and heterozygous animals. These findings were consistent in sensitivity analyses using an alternative distance metric and were not explained by major global differences in multivariate dispersion. Visualisation of the molecular profiles showed relatively stable age-course trajectories in wild-type animals and more dynamic, tissue-dependent changes in heterozygous animals. Gene-specific models provided complementary evidence that individual markers contributed differently to the coordinated multivariate patterns.

Conclusions
This exploratory study shows that small-scale molecular datasets can be analysed as coordinated age-course profiles rather than as independent marker-level outcomes. Integrating constrained permutation-based multivariate testing, visualisation of molecular trajectories and gene-specific follow-up models provides an interpretable framework for studying molecular profiles across tissues and biological conditions. This approach may help bridge classical biostatistical modelling and modern molecular data analysis in experimental studies with limited sample size but complex temporal structure.

Funding: Horizon Europe MSCA-DN programme CHANGE No 101072766.