Load mass spectrometry-based proteomics data using readQFeatures() (original) (raw)
The QFeatures
class
The QFeatures
class stores data as a list of SummarizedExperiment
objects that contain data processed at different levels. For instance, a QFeatures
object may contain data at the peptide-to-spectrum-match (PSM) level, at the peptide level and at the protein level. We call each SummarizedExperiment
object contained in a QFeatures
object a set. Because the different sets are often related, they often share the same samples (columns). QFeatures
automatically creates links between the related samples and their annotations (stored in a singlecolData
table). Similarly, different sets often share related features (rows). For instance, proteins are composed of peptides and peptides are composed of PSMs. QFeatures
automatically creates links between the related features through an AssayLinks
object.
Figure 1: The QFeatures
data class
The QFeatures
object contains a list of SummarizedExperiment
ojects (see class description) on SingleCellExperiment
and QFeatures
objects
library("QFeatures")
Converting tabular data
QFeatures
is designed to process and manipulate the MS-based proteomics data obtained after identification and quantification of the raw MS files. The identification and quantification steps are generally performed by dedicated software (e.g. Sage, FragPipe, Proteome Discoverer, MaxQuant, …) that return a set of tabular data.readQFeatures()
converts these tabular data into a QFeatures
object. We refer to these tables as the assayData
tables.
We distinguish between two use cases: the single-set case and the multi-set case.
The single-set case
The single-set case will generate a QFeatures
object with a singleSummarizedExperiment
object. This is generally the case when reading data at the peptide or protein level, or when the samples where multiplexed (e.g. using TMT) within a single MS run. There are two types of columns:
- Quantitative columns (
quantCols
): 1 to n (depending on technology) - Feature annotations: e.g. peptide sequence, ion charge, protein name
In this case, each quantitative column contains information for a single sample. This can be schematically represented as below:
Figure 2: Schematic representation of a data table under the single-set case
Quantification columns (quantCols
) are represented by different shades of red.
The hyperLOPIT data is an example data that falls under the single-set case (see ?hlpsms
for more details). The quantCols
are X126
,X127N
, X127C
, …, X130N
, X130C
, X131
and correspond to different TMT labels.
In this toy example, there are 3,010 rows corresponding to features (quantified PSMs) and 28 columns corresponding to different data fields generated by MaxQuant during the analysis of the raw MS spectra. The table is converted to a QFeatures
object as follows:
data("hlpsms")
quantCols <- grep("^X", colnames(hlpsms))
(qfSingle <- readQFeatures(hlpsms, quantCols = quantCols))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] quants: SummarizedExperiment with 3010 rows and 10 columns
The object returned by readQFeatures()
is a QFeatures
object containing 1 SummarizedExperiment
set. The set is namedquants
by default, but we could name it psms
by providing thename
argument:
(qfSingle <- readQFeatures(hlpsms, quantCols = quantCols, name = "psms"))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] psms: SummarizedExperiment with 3010 rows and 10 columns
The multi-set case
The multi-set case will generate a QFeatures
object with multipleSummarizedExperiment
objects. This is generally the case when reading data at the PSM level that has been acquired as part of multiple runs. In this case, the identification and quantification software concatenates the results across MS runs in a single table. There are three types of columns:
- Run identifier column (
runCol
): e.g. file name. - Quantification columns (
quantCols
): 1 to n (depending on technology). - Feature annotations: e.g. peptide sequence, ion charge, protein name.
Each quantitative column contains information for multiple samples. This can be schematically represented as below:
Figure 3: Schematic representation of a data table under the multi-set case
Quantification columns (quantCols
) are coloured by run and shaded by label. Every sample is uniquely represented by a colour and shade. Note that every quantCol
contains multiple samples.
We will again use hyperLOPIT data and simulate it was acquired as part of multiple runs, hence falling under the multi-set case. The MS run is often identified with the name of the file it generated.
hlpsms$FileName <- rep(
rep(paste0("run", 1:3, ".raw"), each = 4),
length.out = nrow(hlpsms)
)
Note that the data set now has a column called “FileName” with 3 different runs:
To avoid that a quantification column contains data from multiple samples, readQFeatures()
splits the table into mulitple set depending on the runCol
column, here given as FileName
:
(qfMulti <- readQFeatures(hlpsms, quantCols = quantCols, runCol = "FileName"))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Splitting data in runs.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 3 sets:
#>
#> [1] run1.raw: SummarizedExperiment with 1004 rows and 10 columns
#> [2] run2.raw: SummarizedExperiment with 1004 rows and 10 columns
#> [3] run3.raw: SummarizedExperiment with 1002 rows and 10 columns
The object returned by readQFeatures()
is a QFeatures
object containing 3 SummarizedExperiment
sets. The sets are automatically named based on the values found in runCol
.
Including sample annotations
Data often comes with sample annotations that provide information about the experimental design. These data are generally created by the user. To facilitate sample annotations, readQFeatures()
also allows providing the annotation table as the colData
argument. Depending on the use case, one or multiple columns are required.
For the single-set case, the colData
table must contain a column named quantCols
.
Figure 4: colData
for the single-set case
Let’s simulate such a table:
(coldata <- DataFrame(
quantCols = quantCols,
condition = rep(c("A", "B"), 5),
batch = rep(c("batch1", "batch2"), each = 5)
))
#> DataFrame with 10 rows and 3 columns
#> quantCols condition batch
#> <integer> <character> <character>
#> 1 1 A batch1
#> 2 2 B batch1
#> 3 3 A batch1
#> 4 4 B batch1
#> 5 5 A batch1
#> 6 6 B batch2
#> 7 7 A batch2
#> 8 8 B batch2
#> 9 9 A batch2
#> 10 10 B batch2
We can now provide the table to readQFeatures()
:
(qfSingle <- readQFeatures(hlpsms, quantCols = quantCols, colData = coldata))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] quants: SummarizedExperiment with 3010 rows and 10 columns
For convenience, the quantCols
argument can be omitted when providing colData
(quantCols
are then fetched from this table):
(qfSingle <- readQFeatures(hlpsms, colData = coldata))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] quants: SummarizedExperiment with 3010 rows and 10 columns
The annotations are retrieved as follows:
colData(qfSingle)
#> DataFrame with 10 rows and 3 columns
#> quantCols condition batch
#> <integer> <character> <character>
#> X126 1 A batch1
#> X127C 2 B batch1
#> X127N 3 A batch1
#> X128C 4 B batch1
#> X128N 5 A batch1
#> X129C 6 B batch2
#> X129N 7 A batch2
#> X130C 8 B batch2
#> X130N 9 A batch2
#> X131 10 B batch2
For the multi-set case, the colData
table must contain a column named quantCols
and a column called runCol
.
Figure 5: colData
for the multi-set case
Let’s simulate an annotation table based on our previous example by duplicating the table for each run:
coldataMulti <- DataFrame()
for (run in paste0("run", 1:3, ".raw")) {
coldataMulti <- rbind(coldataMulti, DataFrame(runCol = run, coldata))
}
coldataMulti
#> DataFrame with 30 rows and 4 columns
#> runCol quantCols condition batch
#> <character> <integer> <character> <character>
#> 1 run1.raw 1 A batch1
#> 2 run1.raw 2 B batch1
#> 3 run1.raw 3 A batch1
#> 4 run1.raw 4 B batch1
#> 5 run1.raw 5 A batch1
#> ... ... ... ... ...
#> 26 run3.raw 6 B batch2
#> 27 run3.raw 7 A batch2
#> 28 run3.raw 8 B batch2
#> 29 run3.raw 9 A batch2
#> 30 run3.raw 10 B batch2
We can provide the table to readQFeatures()
:
(qfMulti <- readQFeatures(
hlpsms, quantCols = quantCols, colData = coldataMulti,
runCol = "FileName"
))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Splitting data in runs.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 3 sets:
#>
#> [1] run1.raw: SummarizedExperiment with 1004 rows and 10 columns
#> [2] run2.raw: SummarizedExperiment with 1004 rows and 10 columns
#> [3] run3.raw: SummarizedExperiment with 1002 rows and 10 columns
Additional information
Sample names
readQFeatures()
automatically assigns names that are unique across all samples in all sets. In the single-set case, sample names are provided by quantCols
.
colnames(qfSingle)
#> CharacterList of length 1
#> [["quants"]] X126 X127C X127N X128C X128N X129C X129N X130C X130N X131
In the multi-set case, sample names are the concatenation of the run name and the quantCols (separated by a _
).
colnames(qfMulti)
#> CharacterList of length 3
#> [["run1.raw"]] run1.raw_X126 run1.raw_X127C ... run1.raw_X130N run1.raw_X131
#> [["run2.raw"]] run2.raw_X126 run2.raw_X127C ... run2.raw_X130N run2.raw_X131
#> [["run3.raw"]] run3.raw_X126 run3.raw_X127C ... run3.raw_X130N run3.raw_X131
Special case: empty samples
In some rare cases, it can be beneficial to remove samples where all quantifications are NA
. This can occur when the raw data are searched for labels that were not used during the experiment. For instance, some may quantifying the raw data expecting TMT-16 labelling while the experiment used TMT-11 labels, or used half of the TMT-16 labels. The missing label channels are filled with NA
s. When settingremoveEmptyCols = TRUE
, readQFeatures()
automatically detects and removes columns containing only NA
s.
hlpsms$X126 <- NA
(qfNoEmptyCol <- readQFeatures(
hlpsms, quantCols = quantCols, removeEmptyCols = TRUE
))
#> Checking arguments.
#> Loading data as a 'SummarizedExperiment' object.
#> Formatting sample annotations (colData).
#> Formatting data as a 'QFeatures' object.
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] quants: SummarizedExperiment with 3010 rows and 9 columns
Note that we have set all values in X126
to missing. Hence, the set contains only 9 columns instead of the previous 10.
Reducing verbose
Every call to readQFeatures()
prints progression to the console. To disable the console output, you can use the verbose
argument:
(qfSingle <- readQFeatures(
hlpsms, quantCols = quantCols, verbose = FALSE
))
#> An instance of class QFeatures (type: bulk) with 1 set:
#>
#> [1] quants: SummarizedExperiment with 3010 rows and 10 columns
Under the hood
readQFeatures
proceeds as follows:
- The
assayData
table must be provided as adata.frame
(or any format that can be coerced to adata.frame
).readQFeatures()
converts the table to aSingleCellExperiment
object usingquantCols
to identify the quantitative values that are stored in theassay
slot. Any other column is considered as feature annotation and will be stored asrowData
.
Figure 6: Step1: Convert the input table to a SingleCellExperiment
object
- (Only for the multi-set case:) The
SingleCellExperiment
object is split according to the acquisition run provided by therunCol
column inassayData
.
Figure 7: Step2: Split by acquisition run
- The sample annotations are generated. If no
colData
is provided, the sample annotations are empty. Otherwise,readQFeatures()
matches the information fromassayData
andcolData
based onquantCols
(single-set case) orquantCols
andrunCol
(multi-set case). Sample annotations are stored in thecolData
slot of theQFeatures
object.
Figure 8: Step3: Adding and matching the sample annotations
- Finally, the
SummarizedExperiment
sets and thecolData
are converted to aQFeatures
object.
Figure 9: Step4: Converting to a QFeatures
What about other input formats?
readQFeatures()
should work with any PSM quantification table that is output by a pre-processing software. For instance, you can easily import the PSM tables generated by Proteome Discoverer. The run names are contained in the File ID
column (that should be supplied as therunCol
argument to readQFeatures()
). The quantification columns are contained in the columns starting with Abundance
, eventually followed by a multiplexing tag name. These columns should be stored in a dedicated column in the colData
data to be supplied as runCol
to readQFeatures()
.
The QFeatures
package is meant for both label-free and multiplexed proteomics data. Importing LFQ data is similar to the examples above with the only difference that quantCols
would have only 1 element.
The readSCPfromDIANN()
function is adapted to import label-free and plexDIA/mTRAQ Report.tsv
files generated by DIA-NN.
For more information, see the ?readQFeatures()
and?readQFeaturesFromDIANN()
manual pages, that described the main principle that concern the data import and formatting.
If your input cannot be loaded using the procedure described in this vignette, you can submit a feature request (see next section).
Need help?
You can open an issue on the GitHub repositoryin case of troubles when loading your data with readQFeatures()
. Any suggestion or feature request about the function or the documentation are also warmly welcome.
Session information
R version 4.5.0 (2025-04-11)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.04.2 LTS
Matrix products: default
BLAS: /home/biocbuild/bbs-3.22-bioc/R/lib/libRblas.so
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_GB LC_COLLATE=C
[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: America/New_York
tzcode source: system (glibc)
attached base packages:
[1] grid stats4 stats graphics grDevices utils datasets
[8] methods base
other attached packages:
[1] DT_0.33 ComplexHeatmap_2.25.0
[3] gplots_3.2.0 dplyr_1.1.4
[5] ggplot2_3.5.2 QFeatures_1.19.2
[7] MultiAssayExperiment_1.35.3 SummarizedExperiment_1.39.0
[9] Biobase_2.69.0 GenomicRanges_1.61.0
[11] GenomeInfoDb_1.45.3 IRanges_2.43.0
[13] S4Vectors_0.47.0 BiocGenerics_0.55.0
[15] generics_0.1.4 MatrixGenerics_1.21.0
[17] matrixStats_1.5.0 BiocStyle_2.37.0
loaded via a namespace (and not attached):
[1] tidyselect_1.2.1 farver_2.1.2 bitops_1.0-9
[4] fastmap_1.2.0 lazyeval_0.2.2 digest_0.6.37
[7] lifecycle_1.0.4 cluster_2.1.8.1 Cairo_1.6-2
[10] ProtGenerics_1.41.0 statmod_1.5.0 magrittr_2.0.3
[13] compiler_4.5.0 rlang_1.1.6 sass_0.4.10
[16] tools_4.5.0 igraph_2.1.4 yaml_2.3.10
[19] knitr_1.50 htmlwidgets_1.6.4 S4Arrays_1.9.0
[22] labeling_0.4.3 DelayedArray_0.35.1 plyr_1.8.9
[25] RColorBrewer_1.1-3 abind_1.4-8 KernSmooth_2.23-26
[28] withr_3.0.2 purrr_1.0.4 caTools_1.18.3
[31] colorspace_2.1-1 iterators_1.0.14 scales_1.4.0
[34] gtools_3.9.5 MASS_7.3-65 dichromat_2.0-0.1
[37] tinytex_0.57 cli_3.6.5 rmarkdown_2.29
[40] crayon_1.5.3 rjson_0.2.23 httr_1.4.7
[43] reshape2_1.4.4 BiocBaseUtils_1.11.0 cachem_1.1.0
[46] stringr_1.5.1 parallel_4.5.0 AnnotationFilter_1.33.0
[49] BiocManager_1.30.25 XVector_0.49.0 vctrs_0.6.5
[52] Matrix_1.7-3 jsonlite_2.0.0 bookdown_0.43
[55] GetoptLong_1.0.5 clue_0.3-66 crosstalk_1.2.1
[58] magick_2.8.6 foreach_1.5.2 limma_3.65.1
[61] tidyr_1.3.1 jquerylib_0.1.4 glue_1.8.0
[64] codetools_0.2-20 shape_1.4.6.1 stringi_1.8.7
[67] gtable_0.3.6 UCSC.utils_1.5.0 tibble_3.2.1
[70] pillar_1.10.2 htmltools_0.5.8.1 circlize_0.4.16
[73] R6_2.6.1 doParallel_1.0.17 evaluate_1.0.3
[76] lattice_0.22-7 png_0.1-8 msdata_0.49.0
[79] bslib_0.9.0 Rcpp_1.0.14 SparseArray_1.9.0
[82] xfun_0.52 GlobalOptions_0.1.2 MsCoreUtils_1.21.0
[85] pkgconfig_2.0.3