Overview

LongOmicsTraj separates longitudinal trajectory estimation from trajectory representation and clustering. In this vignette, the packaged GLUCOLD example is used to demonstrate a complete topology-aware clustering workflow:

construct visit-level empirical trajectories;
calculate adjacent changes and assign feature-level Omics Trajectory Signatures (OTS);
cluster genes using combined group-by-visit profiles;
summarise topology composition and purity within clusters;
evaluate topology coverage and duplication;
visualise cluster centres, trajectory arms and topology deviations; and
fit topology-aware FlexMix trajectory clusters in a selected treatment arm.

The example is intentionally small so that it can be run as part of package documentation. The packaged object contains 71 molecular features measured across three visits and four treatment groups.

Required packages

library(LongOmicsTraj)
library(S4Vectors)
library(SummarizedExperiment)
library(MultiAssayExperiment)
library(dplyr)
library(ggplot2)

The optional sections use cluster, matrixStats, patchwork, scales, withr and flexmix.

Load and prepare the GLUCOLD example

data("lot_glucold")

visits <- c(
  "Baseline",
  "6 Months",
  "30 Months"
)

assay_name <- "transcriptomics"

lot <- lot_glucold

experiments <- MultiAssayExperiment::experiments(
  lot@mae
)

se <- experiments[[assay_name]]

metadata <- as.data.frame(
  SummarizedExperiment::colData(se)
)

metadata$visit <- factor(
  metadata$visit,
  levels = visits,
  ordered = TRUE
)

metadata$group <- factor(
  metadata$group
)

stopifnot(
  !anyNA(metadata$visit)
)

SummarizedExperiment::colData(se) <- S4Vectors::DataFrame(
  metadata,
  row.names = rownames(metadata)
)

experiments[[assay_name]] <- se

lot@mae <- MultiAssayExperiment::MultiAssayExperiment(
  experiments = experiments
)

methods::validObject(lot)
#> [1] TRUE

lot
#> 
#> LongOmicsTraj object
#> ────────────────────────────────────────────
#> 
#> Assays
#>   transcriptomics: 71 features × 221 samples
#> 
#> Dimensions
#>   Features : 71
#>   Samples  : 221
#>   Subjects : 90
#> 
#> Study design
#>   Visits : Baseline → 6 Months → 30 Months
#>   Groups : 4
#>            Placebo, ICS 30 months, ICS/LABA 30 months, ICS 6 months then withdrawal
#> 
#> Computed results
#>   trajectories: –
#>   deltas      : –
#>   thresholds  : –
#>   OTS         : –
#> 
#> Estimators
#>   None recorded

The visit order is set explicitly because topology labels depend on the chronological order of adjacent transitions.

Construct feature-level trajectories

For this vignette, empirical group-by-visit means are used. The same clustering workflow can be applied after constructing a VMM with lot_from_lmm(), lot_from_gam() or lot_from_masigpro().

lot_features <- lot_from_empirical(
  object = lot,
  assay = assay_name,
  visits = visits,
  group_col = "group",
  visit_col = "visit",
  estimator = "empirical",
  overwrite = TRUE,
  verbose = FALSE
)

utils::head(
  as.data.frame(lot_features@vmm)
)
#>        object_id object_type           assay   group    visit estimated_value
#> ATP2A3    ATP2A3        gene transcriptomics Placebo Baseline        7.066269
#> CAT          CAT        gene transcriptomics Placebo Baseline        6.852340
#> CCL24      CCL24        gene transcriptomics Placebo Baseline        3.446148
#> CCL26      CCL26        gene transcriptomics Placebo Baseline        3.097151
#> CCL3        CCL3        gene transcriptomics Placebo Baseline        4.701020
#> CCL4        CCL4        gene transcriptomics Placebo Baseline        2.435094
#>        estimator
#> ATP2A3 empirical
#> CAT    empirical
#> CCL24  empirical
#> CCL26  empirical
#> CCL3   empirical
#> CCL4   empirical

Map feature trajectories to OTS labels

Adjacent visit-to-visit changes are calculated once. The example uses the hybrid thresholding strategy only as a reproducible demonstration; the threshold should be selected and examined through sensitivity analysis in substantive applications.

lot_features <- lot_compute_deltas(
  object = lot_features,
  assay = assay_name,
  visits = visits,
  overwrite = TRUE
)

lot_features <- lot_compute_thresholds(
  object = lot_features,
  assay = assay_name,
  method = "hybrid",
  visits = visits,
  overwrite = TRUE,
  verbose = FALSE
)

lot_features <- lot_map_to_ots(
  object = lot_features,
  assay = assay_name,
  visits = visits,
  overwrite = TRUE,
  verbose = FALSE
)

table(
  lot_features@ots$topology_label
)
#> 
#> down_down down_flat   down_up flat_down flat_flat   flat_up   up_down   up_flat 
#>        13        18        36        47        89        32        33        10 
#>     up_up 
#>         6

Each OTS is an ordered pair of transition states because GLUCOLD contains three visits. For example, up_down denotes an increase from baseline to 6 months followed by a decrease from 6 to 30 months.

Cluster genes across all treatment groups

lot_from_clusters() represents each gene by its complete set of group-by-visit means. Clustering therefore reflects combined longitudinal behaviour across all treatment arms rather than a single endpoint or group.

The example uses k-means with four clusters and standardises each gene trajectory before clustering.

lot_clusters <- lot_from_clusters(
  object = lot_features,
  assay = assay_name,
  group_col = "group",
  visit_col = "visit",
  visits = visits,
  k = 4L,
  estimator = "kmeans",
  scale_features = TRUE,
  seed = 1L,
  require_gene_ots = TRUE,
  overwrite = TRUE,
  verbose = FALSE
)

lot_clusters
#> 
#> LongOmicsTraj object
#> ────────────────────────────────────────────
#> 
#> Assays
#>   transcriptomics: 71 features × 221 samples
#> 
#> Dimensions
#>   Features : 71
#>   Samples  : 221
#>   Subjects : 90
#> 
#> Study design
#>   Visits : Baseline → 6 Months → 30 Months
#>   Groups : 4
#>            Placebo, ICS 30 months, ICS/LABA 30 months, ICS 6 months then withdrawal
#> 
#> Computed results
#>   trajectories: ✓
#>   deltas      : ✓
#>   thresholds  : ✓
#>   OTS         : ✓
#> 
#> Estimators
#>   empirical

The cluster memberships, cluster sizes, cluster profiles and clustering diagnostics are stored in the results slot.

as.data.frame(
  lot_clusters@results$cluster_sizes
)
#>   cluster_id n_features
#> 1  cluster_1         34
#> 2  cluster_2         15
#> 3  cluster_3         13
#> 4  cluster_4          9

as.data.frame(
  lot_clusters@results$cluster_diagnostics
)
#>   estimator k n_features n_profiles total_withinss average_silhouette
#> 1    kmeans 4         71         12       301.8066           0.273439
#>   scale_features seed
#> 1           TRUE    1

utils::head(
  as.data.frame(
    lot_clusters@results$cluster_membership
  )
)
#>   child_id parent_id       relationship weight           assay method
#> 1   ATP2A3 cluster_3 cluster_membership      1 transcriptomics kmeans
#> 2      CAT cluster_4 cluster_membership      1 transcriptomics kmeans
#> 3    CCL24 cluster_1 cluster_membership      1 transcriptomics kmeans
#> 4    CCL26 cluster_3 cluster_membership      1 transcriptomics kmeans
#> 5     CCL3 cluster_1 cluster_membership      1 transcriptomics kmeans
#> 6     CCL4 cluster_2 cluster_membership      1 transcriptomics kmeans

The VMM now contains cluster-centre trajectories rather than gene-level trajectories.

utils::head(
  as.data.frame(lot_clusters@vmm)
)
#>   object_id object_type           assay                        group    visit
#> 1 cluster_1     cluster transcriptomics                      Placebo Baseline
#> 2 cluster_1     cluster transcriptomics                ICS 30 months Baseline
#> 3 cluster_1     cluster transcriptomics           ICS/LABA 30 months Baseline
#> 4 cluster_1     cluster transcriptomics ICS 6 months then withdrawal Baseline
#> 5 cluster_1     cluster transcriptomics                      Placebo 6 Months
#> 6 cluster_1     cluster transcriptomics                ICS 30 months 6 Months
#>   estimator estimated_value value_type
#> 1    kmeans      -0.1096570     scaled
#> 2    kmeans      -0.1127952     scaled
#> 3    kmeans       0.2324702     scaled
#> 4    kmeans       1.1709146     scaled
#> 5    kmeans      -0.4987080     scaled
#> 6    kmeans      -0.1756686     scaled

Map cluster centres to topology labels

Cluster centres can be processed by the same delta, threshold and OTS functions used for molecular features.

lot_clusters <- lot_compute_deltas(
  object = lot_clusters,
  assay = assay_name,
  visits = visits,
  overwrite = TRUE
)

lot_clusters <- lot_compute_thresholds(
  object = lot_clusters,
  assay = assay_name,
  method = "hybrid",
  visits = visits,
  overwrite = TRUE,
  verbose = FALSE
)

lot_clusters <- lot_map_to_ots(
  object = lot_clusters,
  assay = assay_name,
  visits = visits,
  overwrite = TRUE,
  verbose = FALSE
)

as.data.frame(
  lot_clusters@ots
)
#>    object_id object_type           assay                        group estimator
#> 1  cluster_1     cluster transcriptomics                ICS 30 months    kmeans
#> 2  cluster_1     cluster transcriptomics ICS 6 months then withdrawal    kmeans
#> 3  cluster_1     cluster transcriptomics           ICS/LABA 30 months    kmeans
#> 4  cluster_1     cluster transcriptomics                      Placebo    kmeans
#> 5  cluster_2     cluster transcriptomics                ICS 30 months    kmeans
#> 6  cluster_2     cluster transcriptomics ICS 6 months then withdrawal    kmeans
#> 7  cluster_2     cluster transcriptomics           ICS/LABA 30 months    kmeans
#> 8  cluster_2     cluster transcriptomics                      Placebo    kmeans
#> 9  cluster_3     cluster transcriptomics                ICS 30 months    kmeans
#> 10 cluster_3     cluster transcriptomics ICS 6 months then withdrawal    kmeans
#> 11 cluster_3     cluster transcriptomics           ICS/LABA 30 months    kmeans
#> 12 cluster_3     cluster transcriptomics                      Placebo    kmeans
#> 13 cluster_4     cluster transcriptomics                ICS 30 months    kmeans
#> 14 cluster_4     cluster transcriptomics ICS 6 months then withdrawal    kmeans
#> 15 cluster_4     cluster transcriptomics           ICS/LABA 30 months    kmeans
#> 16 cluster_4     cluster transcriptomics                      Placebo    kmeans
#>    topology_label n_transitions
#> 1       flat_down             2
#> 2       flat_flat             2
#> 3       flat_flat             2
#> 4       flat_down             2
#> 5         down_up             2
#> 6       flat_flat             2
#> 7         down_up             2
#> 8         flat_up             2
#> 9         up_down             2
#> 10      flat_down             2
#> 11      down_flat             2
#> 12        up_down             2
#> 13        flat_up             2
#> 14      flat_flat             2
#> 15        flat_up             2
#> 16        flat_up             2

These cluster-centre OTS labels describe the representative trajectory of each cluster. They should be interpreted together with the feature-level topology composition stored during clustering.

Inspect topology composition within clusters

When feature-level OTS labels are present before clustering, lot_from_clusters() stores the within-cluster topology distribution in cluster_ots_distribution.

cluster_ots_distribution <- as.data.frame(
  lot_clusters@results$cluster_ots_distribution
)

utils::head(
  cluster_ots_distribution
)
#>   cluster_id                        group topology_label count proportion
#> 1  cluster_1                ICS 30 months      flat_flat    13 0.38235294
#> 2  cluster_1                ICS 30 months      flat_down    12 0.35294118
#> 3  cluster_1                ICS 30 months      down_down     4 0.11764706
#> 4  cluster_1                ICS 30 months        up_down     4 0.11764706
#> 5  cluster_1                ICS 30 months      down_flat     1 0.02941176
#> 6  cluster_1 ICS 6 months then withdrawal      flat_flat    20 0.58823529

Trajectory-arm plot

The trajectory-arm plot shows that a conventional molecular-feature cluster may contain several distinct OTS classes. Arrow width is proportional to the within-cluster topology proportion.

lot_plot_cluster_trajectory_arms(
  object = lot_clusters,
  top_n = 5L,
  min_proportion = 0.05,
  show_counts = TRUE,
  show_percent = TRUE,
  scale_width = TRUE,
  ncol = 2
)

Feature-level trajectory topologies within each cluster and treatment group. Arrow width represents the within-cluster topology proportion.

Quantify topology purity

Topology purity is the proportion of genes assigned to the dominant OTS within each cluster-group combination. Entropy provides a complementary measure of within-cluster topology heterogeneity.

purity_assignments <- as.data.frame(
  lot_clusters@results$gene_ots
) |>
  dplyr::transmute(
    gene_id = as.character(object_id),
    cluster = as.character(cluster_id),
    group = as.character(group),
    topology_label = as.character(topology_label)
  )

purity <- lot_cluster_purity(
  assignments = purity_assignments,
  purity_threshold = 0.90
)

purity$cluster_group_purity
#> # A tibble: 16 × 10
#>    cluster   group       n_genes dominant_topology dominant_count purity entropy
#>    <chr>     <chr>         <int> <chr>                      <int>  <dbl>   <dbl>
#>  1 cluster_1 ICS 30 mon…      34 flat_flat                     13  0.382   1.34 
#>  2 cluster_1 ICS 6 mont…      34 flat_flat                     20  0.588   1.27 
#>  3 cluster_1 ICS/LABA 3…      34 flat_flat                     14  0.412   1.57 
#>  4 cluster_1 Placebo          34 flat_down                     10  0.294   1.80 
#>  5 cluster_2 ICS 30 mon…      15 down_up                        5  0.333   1.45 
#>  6 cluster_2 ICS 6 mont…      15 flat_flat                      7  0.467   1.52 
#>  7 cluster_2 ICS/LABA 3…      15 flat_up                        8  0.533   1.16 
#>  8 cluster_2 Placebo          15 flat_up                        7  0.467   1.62 
#>  9 cluster_3 ICS 30 mon…      13 up_down                        7  0.538   1.16 
#> 10 cluster_3 ICS 6 mont…      13 flat_flat                      4  0.308   1.71 
#> 11 cluster_3 ICS/LABA 3…      13 down_up                        3  0.231   1.80 
#> 12 cluster_3 Placebo          13 flat_down                      8  0.615   0.859
#> 13 cluster_4 ICS 30 mon…       9 flat_up                        4  0.444   1.21 
#> 14 cluster_4 ICS 6 mont…       9 flat_flat                      5  0.556   1.30 
#> 15 cluster_4 ICS/LABA 3…       9 flat_up                        2  0.222   1.89 
#> 16 cluster_4 Placebo           9 down_up                        5  0.556   1.15 
#> # ℹ 3 more variables: n_topologies <int>, normalised_entropy <dbl>,
#> #   passes_purity <lgl>

purity$overall
#> # A tibble: 1 × 6
#>   minimum_purity mean_purity weighted_purity passing_fraction
#>            <dbl>       <dbl>           <dbl>            <dbl>
#> 1          0.222       0.434           0.430                0
#> # ℹ 2 more variables: maximum_normalised_entropy <dbl>,
#> #   n_cluster_group_combinations <int>

A purity of one means that all genes in a cluster-group combination have the same OTS. Lower purity indicates that the cluster combines several trajectory topologies.

Calculate topology coverage and duplication

High purity alone does not guarantee that a clustering solution represents distinct longitudinal behaviours. Several pure clusters may still have the same dominant topology. lot_cluster_topology_score() therefore reports topology coverage and duplicate topology counts separately from purity.

cluster_diagnostics <- as.data.frame(
  lot_clusters@results$cluster_diagnostics
)

minimum_cluster_size <- min(
  as.data.frame(
    lot_clusters@results$cluster_sizes
  )$n_features
)

topology_score <- lot_cluster_topology_score(
  purity = purity,
  k = cluster_diagnostics$k[[1L]],
  minimum_cluster_size = minimum_cluster_size
)

topology_score$group_summary
#> # A tibble: 4 × 13
#>   group               observed_clusters unique_dominant_topo…¹ topology_coverage
#>   <chr>                           <int>                  <int>             <dbl>
#> 1 ICS 30 months                       4                      4              1   
#> 2 ICS 6 months then …                 4                      1              0.25
#> 3 ICS/LABA 30 months                  4                      3              0.75
#> 4 Placebo                             4                      3              0.75
#> # ℹ abbreviated name: ¹unique_dominant_topologies
#> # ℹ 9 more variables: duplicate_topology_count <int>,
#> #   duplicate_topology_fraction <dbl>, minimum_purity <dbl>, mean_purity <dbl>,
#> #   weighted_purity <dbl>, passing_fraction <dbl>,
#> #   maximum_normalised_entropy <dbl>, mean_normalised_entropy <dbl>,
#> #   total_features <int>

topology_score$topology_repetition
#>                           group dominant_topology n_clusters
#> 1                 ICS 30 months           down_up          1
#> 2                 ICS 30 months         flat_flat          1
#> 3                 ICS 30 months           flat_up          1
#> 4                 ICS 30 months           up_down          1
#> 5  ICS 6 months then withdrawal         flat_flat          4
#> 6            ICS/LABA 30 months           flat_up          2
#> 7            ICS/LABA 30 months           down_up          1
#> 8            ICS/LABA 30 months         flat_flat          1
#> 9                       Placebo         flat_down          2
#> 10                      Placebo           down_up          1
#> 11                      Placebo           flat_up          1
#>    is_duplicated_topology
#> 1                   FALSE
#> 2                   FALSE
#> 3                   FALSE
#> 4                   FALSE
#> 5                    TRUE
#> 6                    TRUE
#> 7                   FALSE
#> 8                   FALSE
#> 9                    TRUE
#> 10                  FALSE
#> 11                  FALSE

topology_score$overall
#>   requested_k observed_clusters n_groups n_cluster_group_panels
#> 1           4                 4        4                     16
#>   unique_dominant_topologies unique_group_topology_pairs topology_coverage
#> 1                          5                          11            0.6875
#>   duplicate_topology_count duplicate_topology_fraction minimum_purity
#> 1                        5                      0.3125      0.2222222
#>   mean_purity weighted_purity passing_fraction maximum_normalised_entropy
#> 1   0.4341598       0.4295775                0                  0.9708358
#>   mean_normalised_entropy minimum_cluster_size BIC AIC
#> 1               0.8460981                    9  NA  NA

Topology coverage is the number of unique dominant topologies divided by the number of observed clusters within each group. Values below one indicate that multiple clusters share the same dominant temporal behaviour.

Plot cluster-centre trajectories

Cluster profiles are available in both raw and scaled form.

cluster_profiles <- as.data.frame(
  lot_clusters@results$cluster_profiles
) |>
  dplyr::filter(
    value_type == "scaled"
  )

ggplot(
  cluster_profiles,
  aes(
    x = visit,
    y = value,
    group = group,
    colour = group
  )
) +
  geom_line(
    linewidth = 0.9
  ) +
  geom_point(
    size = 1.8
  ) +
  facet_wrap(
    ~ cluster_id,
    scales = "free_y"
  ) +
  labs(
    title = "Cluster-centre longitudinal profiles",
    x = "Visit",
    y = "Scaled cluster centre",
    colour = "Group"
  ) +
  theme_bw() +
  theme(
    axis.text.x = element_text(
      angle = 30,
      hjust = 1
    ),
    legend.position = "bottom"
  )

Scaled cluster-centre trajectories across treatment groups and visits.

Alternative clustering with PAM

lot_from_clusters() also supports partitioning around medoids. PAM is useful when a medoid-based solution is preferred and the optional cluster package is available.

lot_pam <- lot_from_clusters(
  object = lot_features,
  assay = assay_name,
  group_col = "group",
  visit_col = "visit",
  visits = visits,
  k = 4L,
  estimator = "pam",
  scale_features = TRUE,
  require_gene_ots = TRUE,
  overwrite = TRUE,
  verbose = FALSE
)

as.data.frame(
  lot_pam@results$cluster_diagnostics
)
#>   estimator k n_features n_profiles total_withinss average_silhouette
#> 1       pam 4         71         12             NA          0.2412103
#>   scale_features seed
#> 1           TRUE    1

Topology-aware FlexMix clustering within one treatment arm

lot_from_flexmix_clusters() fits finite mixtures of longitudinal regression trajectories within a selected group. Candidate values of k are compared using BIC, topology purity, the fraction of clusters passing the purity rule and minimum cluster size.

The withdrawal arm is used here because it is expected to contain biologically meaningful non-monotonic responses.

withdrawal_group <- "ICS 6 months then withdrawal"

lot_flexmix <- lot_from_flexmix_clusters(
  object = lot_features,
  assay = assay_name,
  group = withdrawal_group,
  visits = visits,
  k = 2:5,
  group_col = "group",
  visit_col = "visit",
  sample_col = "sample_id",
  trajectory_source = "estimated",
  trajectory_estimator = "empirical",
  model = "quadratic",
  scale_features = TRUE,
  purity_threshold = 0.90,
  required_passing_fraction = 0.80,
  min_cluster_size = 5L,
  nrep = 5L,
  seed = 42L,
  store_fits = FALSE,
  overwrite = TRUE,
  verbose = FALSE
)
#> 2 : * * * * *
#> 3 : * * * * *
#> 4 : * * * * *
#> 5 : * * * * *

as.data.frame(
  lot_flexmix@results$flexmix_model_comparison
)
#>   k trajectory_source trajectory_estimator      BIC      AIC    logLik
#> 2 2         estimated            empirical 472.1373 441.8857 -211.9428
#> 3 3         estimated            empirical 428.2069 381.1488 -176.5744
#> 4 4         estimated            empirical 379.9140 316.0495 -139.0247
#> 5 5         estimated            empirical 336.6366 255.9656 -103.9828
#>   n_parameters minimum_purity mean_purity weighted_purity passing_fraction
#> 2            9      0.4375000   0.4823864       0.5070423                0
#> 3           14      0.4705882   0.5030166       0.5070423                0
#> 4           19      0.4705882   0.5694444       0.5352113                0
#> 5           24      0.4615385   0.5357183       0.5352113                0
#>   minimum_cluster_size passes_topology_rule passes_size_rule acceptable
#> 2                   16                FALSE             TRUE      FALSE
#> 3                   17                FALSE             TRUE      FALSE
#> 4                    9                FALSE             TRUE      FALSE
#> 5                    9                FALSE             TRUE      FALSE

lot_flexmix@results$flexmix_selected_k
#> [1] 5

lot_flexmix@results$flexmix_selection_reason
#> [1] "No candidate satisfied both topology and size criteria; the lowest-BIC solution was selected."

The selected resolution is the smallest candidate satisfying the topology purity and cluster-size criteria. If none satisfies both criteria, the lowest-BIC solution is selected.

Plot the cluster-selection path

lot_plot_cluster_selection(
  object = lot_flexmix,
  metrics = c(
    "BIC",
    "weighted_purity",
    "passing_fraction",
    "minimum_cluster_size"
  ),
  purity_threshold = 0.90,
  required_passing_fraction = 0.80,
  min_cluster_size = 5L,
  show_selected = TRUE
)

Likelihood-based and topology-aware statistics across candidate FlexMix cluster resolutions.

Plot representative cluster trajectories and OTS composition

lot_plot_cluster_topology(
  object = lot_flexmix,
  group = withdrawal_group,
  value_type = "scaled",
  show_features = TRUE,
  max_feature_lines = 30L,
  purity_threshold = 0.90,
  sort_by = "purity",
  seed = 1L
)

Representative FlexMix cluster trajectories and their feature-level OTS composition.

Plot topology deviations

The topology-deviation plot compares the frequency of each OTS within a cluster with its overall frequency among all clustered genes. Positive deviations indicate enrichment and negative deviations indicate depletion.

lot_plot_topology_deviation(
  object = lot_flexmix,
  group = withdrawal_group,
  min_abs_deviation = 0.02,
  top_n = 6L
)

Enrichment or depletion of OTS classes within selected FlexMix clusters relative to the full clustered feature set.

Plot representative genes

Genes assigned to selected FlexMix clusters can be retrieved from the stored assignment table and visualised using the original longitudinal measurements.

assignments <- as.data.frame(
  lot_flexmix@results$flexmix_gene_assignments
)

representative_genes <- assignments |>
  dplyr::filter(
    !is.na(topology_label)
  ) |>
  dplyr::group_by(
    cluster
  ) |>
  dplyr::slice_head(
    n = 1L
  ) |>
  dplyr::ungroup() |>
  dplyr::pull(
    gene_id
  )

representative_genes
#> [1] ATP2A3 CD14   CXCL1  CHI3L1 CAT   
#> 71 Levels: ATP2A3 CAT CCL24 CCL26 CCL3 CCL4 CCL5 CD14 CD68 CD79A CD79B ... ZBTB16

lot_plot_trajectory(
  object = lot_features,
  features = representative_genes,
  assay = assay_name,
  groups = withdrawal_group,
  estimator = "empirical",
  show_subjects = TRUE,
  summary = "mean_se",
  facet_scales = "free_y"
)

Observed longitudinal measurements for one representative gene from each selected FlexMix cluster.

Interpretation

This workflow separates three related but distinct questions:

Do genes have similar numerical profiles?
Addressed by k-means, PAM or FlexMix clustering.
Do members of a cluster share the same temporal topology?
Addressed by topology purity and entropy.
Do different clusters represent distinct temporal behaviours?
Addressed by topology coverage and duplicate topology counts.

A clustering solution may have favourable likelihood or compactness while still merging several OTS classes, or it may split one dominant topology into several clusters. LongOmicsTraj therefore reports conventional clustering diagnostics and topology-aware summaries separately rather than combining them into a single unvalidated score.

Session information

sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 LTS
#> 
#> Matrix products: default
#> BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
#>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
#>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: Etc/UTC
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] stats4    stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#>  [1] BiocStyle_2.41.0            scales_1.4.0               
#>  [3] ggplot2_4.0.3               tibble_3.3.1               
#>  [5] tidyr_1.3.2                 dplyr_1.2.1                
#>  [7] LongOmicsTraj_0.99.0        MultiAssayExperiment_1.39.1
#>  [9] SummarizedExperiment_1.43.0 Biobase_2.73.2             
#> [11] GenomicRanges_1.65.4        Seqinfo_1.3.2              
#> [13] IRanges_2.47.5              S4Vectors_0.51.10          
#> [15] BiocGenerics_0.59.12        generics_0.1.4             
#> [17] MatrixGenerics_1.25.0       matrixStats_1.5.0          
#> [19] rmarkdown_2.32             
#> 
#> loaded via a namespace (and not attached):
#>  [1] utf8_1.2.6          sass_0.4.10         SparseArray_1.13.4 
#>  [4] lattice_0.23-1      digest_0.6.39       magrittr_2.0.5     
#>  [7] RColorBrewer_1.1-3  evaluate_1.0.5      grid_4.6.1         
#> [10] fastmap_1.2.0       jsonlite_2.0.0      Matrix_1.7-6       
#> [13] nnet_7.3-21         BiocManager_1.30.27 purrr_1.2.2        
#> [16] modeltools_0.2-25   jquerylib_0.1.4     abind_1.4-8        
#> [19] cli_3.6.6           rlang_1.3.0         XVector_0.53.0     
#> [22] withr_3.0.3         cachem_1.1.0        DelayedArray_0.39.8
#> [25] yaml_2.3.12         otel_0.2.0          S4Arrays_1.13.2    
#> [28] flexmix_2.3-21      tools_4.6.1         buildtools_1.0.0   
#> [31] vctrs_0.7.3         R6_2.6.1            lifecycle_1.0.5    
#> [34] cluster_2.1.8.3     pkgconfig_2.0.3     gtable_0.3.6       
#> [37] bslib_0.12.0        pillar_1.11.1       data.table_1.18.6.1
#> [40] glue_1.8.1          xfun_0.61           tidyselect_1.2.1   
#> [43] sys_3.4.3           knitr_1.52          farver_2.1.2       
#> [46] patchwork_1.3.2     htmltools_0.5.9     labeling_0.4.3     
#> [49] maketools_1.3.2     compiler_4.6.1      S7_0.2.2

Topology-aware clustering of longitudinal trajectories with LongOmicsTraj