Tutorials › Statistical Methods for Single-Cell RNA Sequencing › Linear Mixture Models for Spot Deconvolution

Spatial Transcriptomics · Tutorial 628 of 700

Linear Mixture Models for Spot Deconvolution

Learn to express reference-based spot deconvolution as a constrained linear model and calculate cell-type mixture estimates from normalized spot counts.

Advanced 10 min read

What You'll Learn

  • Write the linear mixture model using reference profiles and spot-level expression.
  • Explain why library-size scaling and reference calibration matter for proportion estimates.
  • Fit two- and three-cell-type mixtures under nonnegativity and sum-to-one constraints.
  • Calculate a constrained least-squares estimate for noisy measurements.
  • Recognize residual patterns and common interpretation errors.

From Spot Counts to a Mixture Model

A spatial spot can contain RNA from several cell types, while a reference data set provides expression profiles for those types separately. Deconvolution uses the contrast between those reference profiles to estimate how much each type contributes to a spot. The linear mixture model makes this idea explicit: spot expression is represented as a weighted combination of reference profiles.

The weights are meaningful as cell-type proportions only when the measurement scale and references support that interpretation. Raw counts from spots with different sequencing depths are not directly comparable to reference profiles on a different scale. The model below therefore uses expression values after compatible normalization, with the assumptions behind that choice stated explicitly.

Definition: Let \(G\) genes be measured for one spot and let \(K\) reference cell types have calibrated expression profiles on the same scale. Write \(y_g\) for the normalized expression of gene \(g\) in the spot and \(r_{gk}\) for the reference expression of gene \(g\) in cell type \(k\). The linear mixture model is $$ y_g=\sum_{k=1}^{K} r_{gk}p_k+\epsilon_g,\qquad g=1,\ldots,G, $$ where \(p_k\) is the contribution of cell type \(k\) and \(\epsilon_g\) is the discrepancy between the model and the observed expression. In matrix form, \(y=Rp+\epsilon\), where \(R\) is the \(G\times K\) reference-profile matrix. A proportion estimate imposes \(p_k\geq 0\) and \(\sum_{k=1}^{K}p_k=1\).

The symbol \(p_k\) is often called a cell-type proportion, but that interpretation needs care. If reference profiles are calibrated to represent expression per cell and cell types contribute comparable amounts of captured RNA per cell—or if differences in RNA content and capture have been appropriately accounted for—the mixture weights can represent cell-number proportions. If not, the weights more directly describe relative contributions to the measured RNA signal. A model cannot turn RNA contributions into cell-number proportions without the necessary calibration.

For raw counts, library size must also be handled. One convenient formulation is to normalize spot counts to a common scale before fitting and construct reference profiles on that same scale. Alternatively, an exposure factor can be included: the expected raw count for gene \(g\) is proportional to the spot’s sequencing depth times \(\sum_k r_{gk}p_k\). In either formulation, the key requirement is that gene-wise values in the spot and references are comparable. The model is a working approximation: it assumes the selected reference profiles adequately represent the contributing cell types and that a shared set of mixture weights describes the measured genes.

Estimating Proportions Under Constraints

Given \(y\) and \(R\), a standard estimate minimizes the squared discrepancy between the observed expression and the mixture prediction, while restricting the weights to the probability simplex. This means every estimated contribution is nonnegative and all contributions sum to one.

$$ \widehat{p} =\underset{p}{\operatorname{argmin}}\ \|y-Rp\|_2^2 \quad\text{subject to}\quad p_k\geq 0\ \text{for all }k,\qquad \sum_{k=1}^{K}p_k=1. $$

Here \(\|y-Rp\|_2^2=\sum_{g=1}^{G}(y_g-\sum_k r_{gk}p_k)^2\) is the residual sum of squares. If expression measurements have substantially different error variances, a weighted version can give less influence to noisier genes; the weights must be chosen to reflect measurement precision rather than to make a desired cell type appear more abundant. Least squares specifies a fitting criterion, not a complete probability model for counts.

A short derivation clarifies the role of the constraints. For two cell types, set \(p_1=p\) and \(p_2=1-p\). Then the predicted expression is \(r_2+p(r_1-r_2)\), so the unconstrained least-squares estimate is the projection of \(y-r_2\) onto the profile difference \(r_1-r_2\):

$$ \widetilde{p} =\frac{(r_1-r_2)^{\mathsf T}(y-r_2)} {(r_1-r_2)^{\mathsf T}(r_1-r_2)}. $$

This follows by differentiating the squared-error objective \(\|y-r_2-p(r_1-r_2)\|_2^2\) with respect to \(p\), setting the derivative equal to zero, and solving. If the denominator is positive, the constrained estimate is \(\widetilde p\) restricted to the interval from 0 to 1: values below 0 are replaced by 0, and values above 1 by 1. The estimated second proportion is \(1-\widehat p\). With more than two cell types, the same constrained objective is solved jointly; one cannot generally estimate each type separately and then normalize the answers without changing the fit.

Worked Example: Two Cell Types With Two Marker Genes

Suppose two calibrated reference profiles report normalized expression per 10,000 units. For genes 1 and 2, cell type A has profile \(r_A=(80,20)^{\mathsf T}\) and cell type B has profile \(r_B=(20,80)^{\mathsf T}\). A hypothetical spot has \(y=(65,35)^{\mathsf T}\). Assume the spot and references use the same normalization, both profiles are appropriate for this spot, and the measured values are adequately described by a mixture without additional gene-specific shifts.

Let \(p_A=p\), so \(p_B=1-p\), with \(0\leq p\leq 1\). The two gene equations are \(65=80p+20(1-p)=20+60p\) and \(35=20p+80(1-p)=80-60p\). The first gives \(p=(65-20)/60=0.75\); the second gives \(p=(80-35)/60=0.75\). Thus the fitted proportions are 0.75 for A and \(1-0.75=0.25\) for B. Substituting them predicts \(0.75(80)+0.25(20)=65\) and \(0.75(20)+0.25(80)=35\), exactly matching the spot values.

In this hypothetical example, the model attributes 75% of the measured mixture to type A and 25% to type B. Interpreting these as cell-number proportions additionally requires the calibration assumptions described above.

Worked Example: Three Cell Types and the Sum-to-One Constraint

Consider three hypothetical reference profiles for three genes, again on a common normalized scale. Let the rows of \(R\) correspond to genes and the columns to cell types A, B, and C: \(R=\begin{pmatrix}90&10&20\\10&80&20\\10&10&60\end{pmatrix}\). A spot has \(y=(52,33,20)^{\mathsf T}\). Assume compatible normalization and representative, calibrated references. Write the unknown proportions as \(a,b,c\), with \(a+b+c=1\) and each proportion nonnegative.

The first two gene equations are \(52=90a+10b+20c\) and \(33=10a+80b+20c\). Substitute \(c=1-a-b\). The first becomes \(52=20+70a-10b\), or \(7a-b=3.2\). The second becomes \(33=20-10a+60b\), or \(-10a+60b=13\). From the first equation, \(b=7a-3.2\). Substitution into the second gives \(-10a+60(7a-3.2)=13\), so \(410a=205\) and \(a=0.5\). Then \(b=7(0.5)-3.2=0.3\), and \(c=1-0.5-0.3=0.2\).

Check the third gene: \(10(0.5)+10(0.3)+60(0.2)=5+3+12=20\), which matches the observed value. The fitted mixture is therefore 50% A, 30% B, and 20% C. All estimates are nonnegative and sum to one, as required.

Worked Example: Least Squares With Noisy Measurements

Use two hypothetical profiles across three genes: \(r_A=(80,20,40)^{\mathsf T}\), \(r_B=(20,80,40)^{\mathsf T}\), and observed spot expression \(y=(64,35,44)^{\mathsf T}\). Assume a common normalized scale and appropriate references, but allow measurement discrepancies, so the three equations need not have one exact solution. Estimate the mixture by constrained least squares.

Set \(d=r_A-r_B=(60,-60,0)^{\mathsf T}\). With \(p\) the contribution of A, the fitted vector is \(r_B+pd\). The unconstrained projection estimate is \(\widetilde p=d^{\mathsf T}(y-r_B)/(d^{\mathsf T}d)\). Here \(y-r_B=(44,-45,4)^{\mathsf T}\), so the numerator is \(60(44)+(-60)(-45)+0(4)=2640+2700=5340\), and the denominator is \(60^2+(-60)^2=7200\). Therefore \(\widetilde p=5340/7200=0.7417\), which already lies between 0 and 1. The constrained estimates are \(\widehat p_A=0.7417\) and \(\widehat p_B=0.2583\).

The fitted expression is \(r_B+0.7417d\), approximately \((64.50,35.50,40.00)^{\mathsf T}\). The observed-minus-fitted residuals are approximately \((-0.50,-0.50,4.00)^{\mathsf T}\), giving a residual sum of squares of \(0.50^2+0.50^2+4.00^2=16.50\). The third gene has the same reference expression in both cell types, so it does not distinguish their proportions; its discrepancy remains in the residual. This estimate is a least-squares compromise, not an exact reconstruction of every observed value.

Diagnostics, Alternatives, and Limits

A fitted proportion vector is not by itself evidence that the reference mixture describes the spot well. Examine residuals \(y-R\widehat p\) across genes. Large residuals, or groups of residuals that are consistently positive or negative, can indicate an unsuitable reference profile, an unrepresented cell type, normalization mismatch, gene-specific regulation in the tissue, or a failure of the linear approximation. A small overall residual can also hide poor fit for a biologically important subset of genes, so inspect gene-level discrepancies rather than only one summary value.

The constraints are essential to a proportion interpretation. Ordinary unconstrained least squares may return negative weights or weights that do not sum to one. Nonnegative least squares prevents negative weights but, without a sum-to-one constraint or a separate scale parameter, its coefficients need not be proportions. Marker-gene ratio methods can be simple and interpretable when markers are selective and stable, but they generally use less information and can be sensitive to marker noise. Constrained least squares uses several genes jointly, but it does not remove bias caused by poor references or mismatched scales.

Do not interpret a fitted contribution as a causal effect, a direct cell count, or a precise biological truth without supporting calibration and uncertainty analysis. Least-squares estimates are especially limited when a spot has substantial count noise, outlier genes, or expression patterns that change within a cell type. The model also treats the selected reference profiles as fixed inputs; uncertainty in those profiles is not automatically reflected in the fitted proportions. These limitations motivate careful reference construction, scale checks, and residual diagnostics before interpreting spatial patterns.

Common mistake: A learner fits unconstrained least squares and obtains \(p_A=1.12\) and \(p_B=-0.12\), then reports these as “112% A and −12% B.” The arithmetic may be a valid unconstrained fit, but those coefficients are not valid proportions. Enforce nonnegativity and the sum-to-one constraint during fitting; do not present infeasible coefficients as cell fractions.

Check Your Understanding

Answer each question, then compare with the brief solutions below.

  1. For two cell types, why must the second proportion equal one minus the first when the weights represent a complete mixture?
  2. Suppose two profiles for one gene are 30 and 70 normalized units, and the spot value is 50. Under an exact two-type mixture, what are the proportions?
  3. In the three-type example, verify that the fitted proportions predict the second gene value of 33.
  4. In the noisy example, why does the third gene not affect the fitted proportions?
  5. What does a negative coefficient from unconstrained least squares imply about reporting that coefficient as a cell-type proportion?

Solutions. (1) A complete mixture has weights summing to one, so \(p_1+p_2=1\) and \(p_2=1-p_1\). (2) Let \(p\) be the contribution of the type with expression 30. Then \(50=30p+70(1-p)=70-40p\), so \(p=0.5\); the other type also has proportion 0.5. (3) Using \(a=0.5,b=0.3,c=0.2\), the prediction is \(10(0.5)+80(0.3)+20(0.2)=5+24+4=33\). (4) Both types have value 40 for that gene, so its profile difference is zero and the gene provides no information about their relative proportions. (5) It is infeasible as a proportion; refit under nonnegativity and sum-to-one constraints rather than reporting it as a cell fraction.

Clinical Biostats