How to Analyze Your Own Plasma Proteomics Data

Oak
proteomics mass-spectrometry olink data-analysis multi-omics
A huge filter-feeding creature on a misty tidal flat sorting glowing droplets into ranked columns as three luminous rivers merge below it.

By the end of this guide you will have a protein-by-sample matrix with documented normalization, a QC report you can show someone, a set of within-person statistics that separate real change from assay noise, and a join key that lets you line proteins up against your own genome and transcriptome. What you need: the delivered proteomics files (either vendor instrument files plus a FASTA, or a normalized abundance table), roughly 32 GB of RAM and 200 GB of disk if you are reprocessing mass spectrometry raw data, Python 3.11 with pandas/numpy/scipy/statsmodels/scikit-learn, and R with limma and sva. Everything below assumes plasma or serum, which is the sample type nearly every personal profile uses, and everything below is measurement and interpretation. Nothing here is diagnosis. Any protein that looks clinically meaningful belongs in front of a physician who can order a validated assay.

1. Inventory the files and work out which measurement technology you have

Before touching statistics, establish what kind of proteomics you were given, because the two dominant technologies fail in different ways and require different preprocessing. Affinity proteomics (Olink’s proximity extension assay, SomaLogic’s aptamer arrays) delivers a wide, pre-normalized table: for Olink you get a long-format CSV or parquet with columns for SampleID, OlinkID, UniProt, Assay, NPX, LOD, PlateID, and QC_Warning, where NPX is a log2-scale relative unit. Mass spectrometry delivers instrument files (.raw for Thermo, .d directories for Bruker timsTOF), optionally .mzML after conversion, and usually a processed search result such as DIA-NN’s report.tsv plus report.pg_matrix.tsv.

Run a file-level inventory first and record checksums, because you will reprocess at least once and you want to know that the inputs did not change.

find . -type f \( -name '*.raw' -o -name '*.parquet' -o -name '*.tsv' -o -name '*.csv' \) \
  -printf '%p\t%s\n' | sort > manifest.tsv
sha256sum $(cut -f1 manifest.tsv) > checksums.txt

The practical difference: affinity panels measure 1,000 to 11,000 pre-selected proteins with high sensitivity into the low picogram range but give you relative units per assay and inherit the vendor’s normalization decisions. Data-independent acquisition mass spectrometry on neat plasma typically quantifies 400 to 700 protein groups without depletion, in absolute-ish intensity units, dominated by albumin, immunoglobulins and apolipoproteins, and it gives you the peptide-level evidence. We prefer having both, and if you only have one, we prefer the affinity panel for longitudinal tracking because run-to-run technical variance is lower and coverage of low-abundance signaling proteins is far better.

2. Reprocess mass spectrometry raw files yourself

Reprocessing matters because the search settings determine your protein list, and the defaults a vendor used are rarely the ones you want for a single person tracked over time. For DIA data we use DIA-NN, and for data-dependent acquisition we use FragPipe (MSFragger plus IonQuant plus Philosopher). Both are free for non-commercial use and both read Thermo and Bruker files natively, so you can skip msconvert unless you want .mzML for archival.

A library-free DIA-NN run against a reviewed human FASTA looks like this:

diann --f run01.raw --f run02.raw --f run03.raw \
  --fasta UP000005640_9606_reviewed.fasta \
  --fasta-search --predictor --gen-spec-lib \
  --cut 'K*,R*,!*P' --missed-cleavages 1 \
  --var-mods 1 --var-mod 'UniMod:35,15.994915,M' \
  --unimod4 \
  --qvalue 0.01 --matrices --reanalyse \
  --mass-acc 15 --mass-acc-ms1 15 \
  --peak-translation --smart-profiling \
  --relaxed-prot-inf --threads 16 \
  --out report.tsv --out-lib lib.speclib

Three flags carry most of the weight. --reanalyse builds an empirical library from your own runs and requantifies against it, which is match-between-runs and it materially reduces missing values across a longitudinal series. --qvalue 0.01 sets the 1% false discovery rate at precursor level, and you should additionally filter protein groups on PG.Q.Value <= 0.01 yourself rather than trusting the matrix file blindly. --mass-acc should be set from the first pass: run one file with --smart-profiling alone, read the inferred mass accuracy out of the log, then fix it for all subsequent runs so that quantification is comparable across dates.

Use the MaxLFQ-style protein quantities in report.pg_matrix.tsv rather than summing precursor intensities. If you want to check the inference, load report.tsv and count distinct Stripped.Sequence per Protein.Group, because any protein group resting on a single peptide is a candidate for exclusion in a longitudinal comparison. In plasma this removes a surprising number of entries.

3. Build one tidy matrix and one metadata table

Every downstream step gets easier if you commit to a single long-format table plus a sample metadata table, with no protein identifiers hiding in column headers. Keep the matrix in log2 space from this point forward, and never silently mix normalized and raw values.

import pandas as pd, numpy as np

pg = pd.read_csv("report.pg_matrix.tsv", sep="\t")
id_cols = ["Protein.Group","Protein.Ids","Protein.Names","Genes","First.Protein.Description"]
long = pg.melt(id_vars=id_cols, var_name="raw_path", value_name="intensity")
long["sample_id"] = long.raw_path.str.extract(r"(S\d{3}_\d{8})")
long["log2"] = np.log2(long.intensity.replace(0, np.nan))
long = long.dropna(subset=["sample_id"])
long.to_parquet("prot_long.parquet")

The metadata table should carry, at minimum, collection date and time of day, fasting duration, plate or batch identifier, instrument and column change dates, and anything you know about the two days before the draw (illness, hard exercise, alcohol, travel, vaccination). Time of day and recent exercise move real plasma proteins by amounts comparable to the effects people go looking for, so an unlabeled sample is close to uninterpretable in a within-person design. For Olink data, carry QC_Warning and LOD through as columns rather than filtering early, so you can test how sensitive a conclusion is to the floor.

4. Run QC before any statistics

QC exists to catch the failure modes that produce confident, wrong answers later. Compute five things for every sample and look at all five, because each catches a different problem.

d = pd.read_parquet("prot_long.parquet")
wide = d.pivot_table(index="Genes", columns="sample_id", values="log2")

qc = pd.DataFrame({
    "n_quant": wide.notna().sum(),
    "median_log2": wide.median(),
    "iqr_log2": wide.quantile(.75) - wide.quantile(.25),
})
corr = wide.dropna().corr(method="spearman")
qc["median_corr"] = corr.median()
print(qc.sort_values("n_quant"))

Identification count per run should be stable within about 10% across a series. A single run 30% low usually means an injection or chromatography problem, not biology. Median correlation against the other samples, computed on complete cases, should sit above roughly 0.95 for technical replicates of the same plasma and above 0.9 for different draws from the same person. Then check the contamination markers explicitly, because plasma is easy to spoil at collection: hemoglobin subunits (HBB, HBA1) indicate hemolysis, and platelet proteins (PF4, PPBP, THBS1) indicate poor centrifugation or delayed processing. A sample with hemoglobin two or three log2 units above your series median should be flagged and probably dropped, since hemolysis inflates hundreds of intracellular proteins at once and will look like a dramatic systemic event.

Finally, plot missingness against mean intensity. In mass spectrometry the two should be strongly coupled, with missing values concentrated in low-abundance proteins. If missingness is instead concentrated in one run date, you have a run problem. If it is random with respect to intensity, check your search settings.

5. Normalize, then decide about imputation on purpose

Normalization removes the part of the between-sample variation that comes from loading, injection volume, and instrument response. Median centering in log2 space handles most of it for plasma and is hard to get wrong. Quantile normalization is stronger and appropriate when intensity distributions differ in shape, at the cost of forcing every sample to the same distribution, which can suppress a genuine global shift. Statistical toolchains for quantitative omics have converged on this general sequence of log transform, normalize, filter, then model, and the reasoning behind the ordering is laid out in the original DAnTE description.1

med = wide.median(axis=0)
norm = wide.sub(med - med.median(), axis=1)
keep = norm.notna().mean(axis=1) >= 0.7          # present in >=70% of samples
norm = norm.loc[keep]

For Olink data, skip this step. NPX is already normalized per plate and per assay, and re-centering it across assays is meaningless because each assay has its own arbitrary intercept. What you may need is bridging normalization if samples were run across plates or across product versions, using the vendor’s bridging samples and a per-assay median offset.

On imputation we are conservative. Filter to proteins quantified in at least 70% of samples and model the remaining gaps rather than filling them, because most mass spectrometry missingness is left-censored (missing not at random) and single-value imputation with a low constant manufactures differences whose size depends on how many values you filled. If a tool requires a complete matrix, use a downshifted-normal draw (mean shifted 1.8 standard deviations down, width 0.3, per sample) and record that the resulting p-values for sparse proteins are not trustworthy. Better: use limma on the observed values, which handles unbalanced missingness without invention.

6. Remove batch structure without removing your own biology

Batch effects in personal proteomics are almost always confounded with time, because samples collected in 2026 were run in 2026. This is the single most dangerous structure in a longitudinal dataset, and no algorithm resolves it. The fix is design: keep aliquots frozen at −80 °C and run them together, or include a shared reference aliquot on every plate or in every run block so that you can estimate the block offset from a sample whose biology is constant.

If you have reference samples, estimate and subtract their per-protein block offsets. If you have only a batch label that is not perfectly confounded with the question you care about, use limma for visualization and sva::ComBat with caution for downstream modeling.

library(limma)
mat <- as.matrix(read.csv("norm_matrix.csv", row.names = 1))
meta <- read.csv("metadata.csv")
vis <- removeBatchEffect(mat, batch = meta$plate,
                         design = model.matrix(~ meta$condition))

Use removeBatchEffect output for plots only. For inference, put the batch term in the model instead, so that the uncertainty it absorbs is carried into the standard errors rather than discarded.

7. Do statistics that make sense for one person

Most proteomics statistics assume two groups of many people. You have many timepoints of one person, so the question changes from “does this protein differ between groups” to “is this value unusual for me, and is this trajectory real”. Three analyses cover most of what you will want.

First, a within-person reference interval. For each protein, compute your own median and a robust scale estimate across baseline timepoints, then express every new sample as a robust z score against your own history. Use median absolute deviation scaled by 1.4826 rather than standard deviation, so one bad draw does not widen the interval.

base = norm[baseline_samples]
mu  = base.median(axis=1)
sd  = 1.4826 * (base.sub(mu, axis=0)).abs().median(axis=1)
z   = norm.sub(mu, axis=0).div(sd.replace(0, np.nan), axis=0)

Second, a trend test. Fit an ordinary least squares or robust regression of log2 abundance on time per protein, then adjust the p-values with Benjamini-Hochberg across proteins. With 500 to 3,000 proteins and 6 to 20 timepoints, expect low power for anything but strong monotone trends, and expect several hits driven entirely by one influential point. Always plot the top 20 trajectories before believing any of them.

Third, a variance decomposition. Run duplicate aliquots of at least two draws so you can estimate technical variance per protein, then compare within-person biological variance against it. Proteins whose biological variance does not exceed technical variance are not trackable for you at this assay’s precision, and knowing which ones those are saves you from chasing noise for years.

Resist the urge to train a classifier on a single person’s timecourse. Machine learning on proteomic data overfits readily when samples are few and correlated, and the practical guardrails (nested cross-validation, held-out data, reporting feature stability rather than a single ranked list) are spelled out clearly in the OmicLearn work.2 With n = 1 there is no defensible held-out set, so descriptive statistics with explicit uncertainty are the right tool.

8. Annotate and interpret, carefully

Interpretation begins with getting identifiers right. Map every protein group to a UniProt accession and a HGNC gene symbol, keep the accession as the primary key, and store the full Protein.Ids string so that ambiguous groups stay visible. Never join on gene symbol alone, since symbols change and several map many-to-one.

For pathway-level reading, use over-representation or rank-based enrichment against Reactome and Gene Ontology, and treat the result as a hypothesis-organizing summary rather than a finding. Enrichment on plasma proteins carries a specific trap: the plasma proteome is a mixture of proteins secreted on purpose (complement, coagulation, apolipoproteins) and proteins that leaked from tissue, so an enriched “cytoskeleton” term often means cell lysis rather than a cytoskeletal program. Methods for guiding pathway and network interpretation with more than one omic at a time, which reduce this ambiguity by requiring concordance, are reviewed in recent protocol work on clinical metabolomics and proteomics.3

9. Join proteins to your genome and transcriptome

The reason to hold proteomics inside a larger personal profile is that each layer constrains the interpretation of the others, and general strategies for that integration are well described.4 Three joins are worth doing immediately.

Start with cis-pQTLs. Pull your genotypes at known protein quantitative trait loci from your VCF and check whether an unusual protein level has a simple genetic explanation, which is common for affinity assays where a coding variant can alter epitope binding and shift the measured value without changing protein abundance.

bcftools query -f '%CHROM\t%POS\t%ID\t%REF\t%ALT[\t%GT]\n' \
  -R pqtl_targets.bed sample.vcf.gz > pqtl_genotypes.tsv

Next, correlate protein against transcript for the genes present in both layers. Expect modest agreement: mRNA and protein correlate imperfectly in general, and for plasma proteins measured against whole-blood RNA the tissue of origin often differs entirely, so a mismatch is informative rather than an error. Finally, align proteins with your continuous glucose and blood biomarker series on a common time axis, since a protein shift that coincides with a documented infection or a training block is interpretable in a way that an isolated value is not. The broader computational problem of reconciling layers measured on different scales and with different noise models is an old one in this field and worth reading about before building your own integration.5

Common problems

Protein groups are not genes. DIA-NN and MSFragger report groups of indistinguishable proteins, and collapsing them to a single symbol silently merges isoforms and family members. Keep the group as the unit of analysis and report it.

Missingness drives false positives. A protein detected in five later samples and none of the earlier ones will show a huge “increase” under any imputation scheme. Check the raw detection pattern for every top hit before writing it down.

Hemolysis masquerades as systemic change. Screen HBB and HBA1 in every sample, every time. One mishandled tube can dominate a principal component.

Dynamic range beats depth. Albumin and immunoglobulins are roughly half of plasma protein mass, so undepleted mass spectrometry will keep returning the same few hundred abundant proteins no matter how long you run the gradient. If you need low-abundance signaling proteins, the answer is an affinity panel or depletion, not more instrument time.

Assay floors are real. Olink values below LOD are not quantitative, and a “change” from 0.3 to 0.8 NPX in a protein whose LOD sits at 1.2 is measurement noise. Track what fraction of your timepoints for a given assay sit above LOD before interpreting any trend.

Batch is confounded with time. If each draw was processed on the day it was collected, you cannot separate assay drift from biology. Bank aliquots and run them in blocks with shared references.

Tool sprawl wastes weeks. The number of published omics tools far exceeds the number that are maintained, and directories built to index them illustrate the scale of the problem.6 Pick one search engine, one statistical framework, and one plotting layer, pin the versions in a lockfile, and record the exact command lines in the repository next to the results.

Anything that looks clinically relevant needs a clinical assay. Research-grade proteomics is not a diagnostic test, units are relative, and reference intervals are your own rather than population-validated. Bring the observation and the raw data to a physician, who can order a validated measurement of the same analyte.

Oak builds longitudinal molecular profiles of individuals: whole-genome sequencing, RNA sequencing, proteomics, blood biomarkers, and continuous glucose data, integrated into one model of you. Build your profile.

Footnotes

  1. Ashoka D. Polpitiya, Wei-Jun Qian, Navdeep Jaitly, et al. DAnTE: a statistical tool for quantitative analysis of -omics data. Bioinformatics, 2008. https://doi.org/10.1093/bioinformatics/btn217 ↩

  2. Furkan M. Torun, Sebastian Virreira Winter, Sophia Doll, et al. Transparent Exploration of Machine Learning for Biomarker Discovery from Proteomics and Omics Data. Journal of Proteome Research, 2022. https://doi.org/10.1021/acs.jproteome.2c00473 ↩

  3. Christina Schmidt, Thomas Naake. Multi-Omics Guided Pathway and Network Analysis of Clinical Metabolomics and Proteomics Data. Methods in Molecular Biology, 2026. https://doi.org/10.1007/978-1-0716-5452-1_17 ↩

  4. Efi Athieniti, George M. Spyrou. A guide to multi-omics data collection and integration for translational medicine. Computational and Structural Biotechnology Journal, 2022. https://doi.org/10.1016/j.csbj.2022.11.050 ↩

  5. Bonnie Berger, Jian Peng, Mona Singh. Computational solutions for omics data. Nature Reviews Genetics, 2013. https://doi.org/10.1038/nrg3433 ↩

  6. V. J. Henry, A. E. Bandrowski, A.-S. Pepin, et al. OMICtools: an informative directory for multi-omic data analysis. Database, 2014. https://doi.org/10.1093/database/bau069 ↩