Import data

The importCustomAtlas() function is the universal entry point. It accepts STAR ReadsPerGene.out.tab files from any plant species:

library(OrnTx)

# Import Rosa chinensis expression data
se <- importCustomAtlas(
  counts_dir    = "path/to/star/counts/",
  metadata_file = "path/to/SraRunTable.txt",
  species       = "Rosa chinensis",
  strandedness  = "reverse",
  min_samples   = 10
)
se

For demonstration, we use the built-in example dataset:

# Load built-in example data
example_files <- system.file("extdata", package = "OrnTx")
cat("OrnTx installed at:", example_files, "\n")

## OrnTx installed at:

# Show built-in annotation
ann <- getRosaAnnotation()
cat("Rosa chinensis genes annotated:", nrow(ann), "\n")

## Rosa chinensis genes annotated: 50455

head(ann)

##     gene_id               locus_name   gene_biotype
## 1 gene49924 RchiOBHm_Chr0c01g0497671 protein_coding
## 2 gene49925 RchiOBHm_Chr0c01g0497681 protein_coding
## 3 gene49926 RchiOBHm_Chr0c01g0497691 protein_coding
## 4 gene49927 RchiOBHm_Chr0c01g0497701           tRNA
## 5 gene49928 RchiOBHm_Chr0c01g0497711 protein_coding
## 6 gene49929 RchiOBHm_Chr0c01g0497721          ncRNA

Normalize expression data

data("rosa_example", package = "OrnTx")
se <- rosa_example
# CPM normalization
se_norm <- normalizeCounts(se, method = "CPM")

## [OrnTx] Normalization complete: assay 'CPM' added.

se_norm

## class: SummarizedExperiment 
## dim: 500 6 
## metadata(6): species strandedness ... import_date ornAtlas_version
## assays(2): counts CPM
## rownames(500): gene10715 gene21520 ... gene37167 gene26373
## rowData names(1): gene_id
## colnames(6): SRR10020018 SRR10020021 ... SRR10604033 SRR10604034
## colData names(126): Run Assay.Type ... tissue_clean tissue_simple

PCA visualization

# PCA of all 270 samples colored by tissue type
p <- plotPCA(se_norm,
             assay_name   = "CPM",
             color_by     = "tissue_simple",
             label_points = FALSE)
print(p)

The PCA reveals distinct tissue-specific transcriptional programs:

• PC1 (38.1%) — captures petal developmental variation
• PC2 (12.3%) — separates floral organ types
• Six tissue types cluster distinctly: petal, abscission zone, bud, meristem, stamen, pistil

Build expression atlas

# Aggregate replicates into condition-level atlas
atlas <- buildAtlas(se_norm,
                    group_var  = "tissue_simple",
                    assay_name = "CPM")

## Building atlas for groups: petal, abscission zone

## Genes retained after filtering: 500

# Heatmap of top 50 variable genes
plotAtlas(atlas,
          top_n  = 50,
          title  = "Rosa chinensis Expression Atlas")

Differential expression analysis

# DE analysis: petal vs abscission zone
results <- runDE(
  se,
  design        = ~ tissue_simple,
  contrast      = c("tissue_simple", "petal", "abscission zone"),
  alpha         = 0.05,
  lfc_threshold = 1
)

## [OrnTx] Running DESeq2 | design: ~tissue_simple

##   Note: levels of factors in the design contain characters other than
##   letters, numbers, '_' and '.'. It is recommended (but not required) to use
##   only letters, numbers, and delimiters '_' or '.', as these are safe characters
##   for column names in R. [This is a message, not a warning or an error]

## [OrnTx] Genes after filtering: 494

##   Note: levels of factors in the design contain characters other than
##   letters, numbers, '_' and '.'. It is recommended (but not required) to use
##   only letters, numbers, and delimiters '_' or '.', as these are safe characters
##   for column names in R. [This is a message, not a warning or an error]

## [OrnTx] Up: 93 | Down: 35 | Total sig: 128

# Summary
cat("Upregulated in petal:", sum(results$direction == "Up"), "\n")

## Upregulated in petal: 93

cat("Downregulated in petal:", sum(results$direction == "Down"), "\n")

## Downregulated in petal: 35

# Volcano plot
plotVolcano(results,
            title = "Petal vs Abscission Zone - Rosa chinensis",
            top_n = 20)

Trait Annotation

library(OrnTx)
data("rosa_example", package = "OrnTx")
results <- runDE(rosa_example,
                 design   = ~ tissue_simple,
                 contrast = c("tissue_simple", "petal", "abscission zone"))

## [OrnTx] Running DESeq2 | design: ~tissue_simple

##   Note: levels of factors in the design contain characters other than
##   letters, numbers, '_' and '.'. It is recommended (but not required) to use
##   only letters, numbers, and delimiters '_' or '.', as these are safe characters
##   for column names in R. [This is a message, not a warning or an error]

## [OrnTx] Genes after filtering: 494

##   Note: levels of factors in the design contain characters other than
##   letters, numbers, '_' and '.'. It is recommended (but not required) to use
##   only letters, numbers, and delimiters '_' or '.', as these are safe characters
##   for column names in R. [This is a message, not a warning or an error]

## [OrnTx] Up: 93 | Down: 35 | Total sig: 128

results_ann <- annotateWithTraits(results, trait = "all")

## [OrnTx] Trait genes retrieved: 64

## [OrnTx] Trait genes: 0

## [OrnTx] Significant trait genes: 0

# Show top significant DE genes
top_genes <- results_ann[results_ann$significant == TRUE, ]
top_genes <- top_genes[order(top_genes$padj), ]
head(top_genes[, c("gene_id", "log2FoldChange", "padj", "direction")], 10)

##             gene_id log2FoldChange         padj direction
## gene33157 gene33157       7.679097 5.661335e-15        Up
## gene13141 gene13141       6.353149 5.355207e-11        Up
## gene44690 gene44690       6.113924 3.554772e-10        Up
## gene25897 gene25897      -4.641978 7.922341e-09      Down
## gene37347 gene37347       6.006105 8.356386e-09        Up
## gene15276 gene15276       7.064469 1.000981e-08        Up
## gene50496 gene50496      10.802829 1.903646e-08        Up
## gene17847 gene17847       4.024180 1.132216e-07        Up
## gene29604 gene29604      11.353449 1.132216e-07        Up
## gene4585   gene4585       9.192342 1.665915e-07        Up

Trait gene analysis

OrnTx includes a curated database of 64 floriculture trait genes anchored to Arabidopsis thaliana orthologs, enabling cross-species trait queries regardless of genome annotation quality.

# Get all flower color genes
color_genes <- getFlowerGenes(trait = 'flower_color')

## [OrnTx] Trait genes retrieved: 17

nrow(color_genes)

## [1] 17

# Get scent genes with strong evidence in Rosa chinensis
rose_scent <- getFlowerGenes(
  trait    = 'scent',
  species  = 'Rosa chinensis',
  evidence = 'Strong'
)

## [OrnTx] Trait genes retrieved: 9

# Find orthologs for known genes
findOrthologs(c('CHS', 'DFR', 'ANS', 'PAL1', 'ACO1'),
              id_type = 'symbol')

## [OrnTx] Queried: 5 | Matched: 5

##    gene_symbol arabidopsis_id                    trait           trait_category
## 1          CHS      AT5G13930 flower_color_anthocyanin             Flower Color
## 5          DFR      AT5G42800 flower_color_anthocyanin             Flower Color
## 6          ANS      AT4G22880 flower_color_anthocyanin             Flower Color
## 18        PAL1      AT3G53260          scent_benzenoid             Floral Scent
## 31        ACO1      AT2G19590      senescence_ethylene Senescence and Vase Life
##                                pathway rose_evidence marigold_evidence key_pmid
## 1             Anthocyanin biosynthesis        Strong          Moderate 15923344
## 5             Anthocyanin biosynthesis        Strong          Moderate 15923344
## 6             Anthocyanin biosynthesis        Strong          Moderate 15923344
## 18 Benzenoid/Phenylpropanoid volatiles        Strong              Weak 17353359
## 31 Ethylene biosynthesis and signaling        Strong          Moderate 20952626

# Standardize developmental stage names across studies
raw_stages <- c('bud', 'Bud', 'floral bud', 'Stage 1',
                'open flower', 'OF', 'full bloom', 'senescent')
standardizeStage(raw_stages)

## [OrnTx] Standardized 6 of 8 stage names

## [1] "bud"         "bud"         "bud"         "bud"         "open_flower"
## [6] "open_flower" "open_flower" "senescent"

Example: Tagetes erecta (Marigold)

Marigold is India’s most economically important floriculture crop, widely used in religious ceremonies, festivals, and as a source of lutein for pharmaceutical production.

# Import marigold expression data
se_marigold <- importCustomAtlas(
  counts_dir    = "path/to/marigold/counts/",
  metadata_file = "path/to/marigold_metadata.csv",
  species       = "Tagetes erecta",
  sample_col    = "Run",
  strandedness  = "reverse",
  min_samples   = 3
)

# Normalize
se_mg_norm <- normalizeCounts(se_marigold, method = "CPM")

# PCA colored by developmental stage and cultivar
p_mg <- plotPCA(se_mg_norm,
                assay_name   = "CPM",
                color_by     = "dev_stage",
                shape_by     = "flower_color",
                label_points = FALSE)
print(p_mg)

# DE: orange (V-01) vs yellow mutant (V-01M) at Stage IV
results_mg <- runDE(
  se_marigold,
  design        = ~ flower_color,
  contrast      = c("flower_color", "orange", "yellow"),
  alpha         = 0.05,
  lfc_threshold = 1
)

plotVolcano(results_mg,
            title = "Orange vs Yellow - Tagetes erecta Stage IV",
            top_n = 15)

Required input format

OrnTx accepts output from the STAR aligner with --quantMode GeneCounts. Each sample generates a ReadsPerGene.out.tab file:

N_unmapped      1234    1234    1234
N_multimapping  5678    5678    5678
N_noFeature     9012    9012    9012
N_ambiguous     3456    3456    3456
gene1           0       0       25
gene2           15      8       7
gene3           234     118     116

The metadata CSV must have a column of sample IDs matching the count file names, plus any experimental variables:

Run,species,tissue,dev_stage,cultivar,replicate
SRR001,Rosa chinensis,petal,Stage_I,Old_Blush,1
SRR002,Rosa chinensis,petal,Stage_I,Old_Blush,2
SRR003,Rosa chinensis,leaf,mature,Old_Blush,1