Bioensayos — Ecuaciones de MatLens

Esta sección describe lo que calcula el módulo de bioensayos. Los cálculos están en matlens/core/bioassays.py; los diálogos de matlens/ui/bioassays_ui.py solo muestran sus resultados. Cuando no hay hoja activa, los diálogos usan un conjunto de datos de ejemplo fijo y lo indican en el título del gráfico.

Ajuste no lineal común

Código: matlens/core/bioassays.py::fit_nls

Mínimos cuadrados no ponderados (scipy.optimize.least_squares, Levenberg–Marquardt o región de confianza con límites):

\min_{\boldsymbol\theta}\ \mathrm{RSS} = \sum_i [y_i - f(x_i;\boldsymbol\theta)]^2,\qquad s^2 = \frac{\mathrm{RSS}}{n-k},\qquad \mathbf C = s^2(\mathbf J^\top\mathbf J)^{-1}

con $\mathbf J$ por diferencias centrales en el óptimo. Intervalos: $\theta \pm t_{1-\alpha/2,\,n-k}\,u(\theta)$ , con el nivel de confianza pedido. $R^2 = 1 - \mathrm{RSS}/\mathrm{TSS}$ . Funciones de los parámetros (EC50 del 5PL, tiempo de duplicación, Cq por SDM): método delta, $u^2(g) = \nabla g^\top\mathbf C\,\nabla g$ .

\mathrm{AICc} = n\ln\frac{\mathrm{RSS}}{n} + 2K + \frac{2K(K+1)}{n-K-1},\quad K = k+1,\qquad w_i = \frac{e^{-\Delta_i/2}}{\sum_j e^{-\Delta_j/2}}

Modelos anidados: prueba F de suma de cuadrados extra, $F = \dfrac{(\mathrm{RSS}_0-\mathrm{RSS}_1)/(\nu_0-\nu_1)}{\mathrm{RSS}_1/\nu_1}$ .

Validación: reproduce los valores certificados NIST StRD Rat42 (4PL) y Rat43 (5PL), incluidos los errores estándar.

Referencias: Burnham y Anderson, Model Selection and Multimodel Inference, 2.ª ed. (2002); Motulsky y Christopoulos, Fitting Models to Biological Data (OUP 2004).

Placas: Z′, SSMD, %CV y atípicos

Código: matlens/core/bioassays.py::plate_qc, z_prime, ssmd, cv_percent, modified_z

Z' = 1 - \frac{3\,(s_p + s_n)}{\lvert \bar x_p - \bar x_n \rvert},\qquad \beta_\text{SSMD} = \frac{\bar x_p - \bar x_n}{\sqrt{s_p^2 + s_n^2}},\qquad \%CV = 100\,\frac{s}{\bar x}

M_i = 0.6745\,\frac{x_i - \tilde x}{\mathrm{MAD}},\qquad \text{atípico si } |M_i| > 3.5

Todas las desviaciones estándar son muestrales ( $n-1$ ). Z′ y SSMD exigen al menos dos pocillos de cada control; si faltan, no se muestra ningún valor y se explica el motivo. Los atípicos se buscan dentro de cada grupo de rol con al menos 3 pocillos ( $\tilde x$ mediana, MAD mediana de las desviaciones absolutas). Clasificación de Z′: $\ge 0.5$ excelente; $0 < Z' < 0.5$ aceptable; $\le 0$ no apto.

Referencias: Zhang, Chung y Oldenburg, J. Biomol. Screen. 4 (1999) 67; Zhang, Genomics 89 (2007) 552; Iglewicz y Hoaglin, How to Detect and Handle Outliers (ASQC 1993).

Viabilidad celular

Código: matlens/core/bioassays.py::viability_percent

\%V = 100\,\frac{S - \bar B}{\bar V - \bar B},\qquad u^2(\%V) = \Big(\frac{\partial\%V}{\partial\bar V}\Big)^2 u^2(\bar V) + \Big(\frac{\partial\%V}{\partial\bar B}\Big)^2 u^2(\bar B)

$\bar B$ es la media de los blancos (0 %), $\bar V$ la de los pocillos de vehículo o, si no hay, del control (−) sin tratar (100 %); $u$ son los errores estándar de las medias (propagación de primer orden, GUM). La señal de cada pocillo se toma como exacta.

Referencia: Riss et al., Cell Viability Assays, Assay Guidance Manual (NCBI, 2013/2016).

Western blot: normalización a la proteína de carga

Código: matlens/core/bioassays.py::western_normalize

r_j = \frac{T_j}{L_j},\qquad \text{relativo}_j = \frac{r_j}{\overline{r}_\text{control}}

$T$ volumen de la banda diana, $L$ de la proteína de carga (o proteína total). Se informa también el %CV de los carriles control.

Referencia: Taylor y Posch, BioMed Res. Int. 2014, 361590.

Dosis-respuesta 4PL / 5PL (IC50 / EC50)

Código: matlens/core/bioassays.py::fit_dose_response, logistic4, logistic5

Con $x = \log_{10} C$ ( $C>0$ ):

y_\text{4PL} = B + \frac{T - B}{1 + 10^{(\log_{10} C_{50} - x)\,h}},\qquad y_\text{5PL} = B + \frac{T - B}{\left[1 + 10^{(\log_{10} C_{50} - x)\,h}\right]^{s}}

El EC50 informado es la concentración de respuesta media entre $B$ y $T$ : en el 4PL, $C_{50}$ ; en el 5PL, $\log_{10}\mathrm{EC}_{50} = \log_{10}C_{50} - \log_{10}(2^{1/s}-1)/h$ . Intervalo por el método delta en $\log_{10}\mathrm{EC}_{50}$ con $t_{n-k}$ y transformado a concentración (asimétrico):

\mathrm{EC}_{50}^{\pm} = 10^{\log_{10}\mathrm{EC}_{50} \pm t\,u(\log_{10}\mathrm{EC}_{50})}

Opciones: $B=0$ o $T=100$ fijos; normalización previa a 0–100 % (desactivada por defecto). Valores iniciales: $B=\min y$ , $T=\max y$ , $\log_{10}C_{50}$ = mediana de $x$ , $h = \pm 4/\mathrm{rango}(x)$ según el sentido de la curva, $s=1$ (acotado a $[10^{-3}, 10^3]$ ).

Comparación de curvas (compare_curves): ajuste global con parámetros comunes (pendiente de Hill: paralelismo; $\log_{10}$ EC50: potencia; o todos) frente a ajustes independientes, con la prueba F de suma de cuadrados extra y el AICc de ambos.

Referencias: Gottschalk y Dunn, Anal. Biochem. 343 (2005) 54; Motulsky y Christopoulos (2004).

Ecuación de Hill

Código: matlens/core/bioassays.py::fit_hill

y = y_0 + y_\text{max}\,\frac{x^n}{K^n + x^n}

con $y_0 = 0$ salvo que se pida línea base. $n$ y $K$ con intervalos $t$ . Es el 4PL con $B = y_0$ , $h = n$ y $\mathrm{EC}_{50} = K$ .

Referencia: Hill, J. Physiol. 40 (1910) iv.

Curva estándar (ELISA / BCA / Bradford)

Código: matlens/core/bioassays.py::fit_standard_curve, back_calculate

Modelos candidatos: lineal $a + bx$ , cuadrático $a + bx + cx^2$ , 4PL $d + (a-d)/[1+(x/c)^b]$ y 5PL $d + (a-d)/[1+(x/c)^b]^g$ (en concentración, admite el estándar cero). Se ajustan todos; se muestra para cada uno $k$ , DE residual, $R^2$ , AICc y peso de Akaike, y se usa el de menor AICc salvo que el usuario elija otro. Si la prueba de rachas de los residuos del modelo usado da $p < 0.05$ se avisa de un patrón sistemático.

Retrocálculo: $x_0$ con $f(x_0) = y_0$ dentro del intervalo calibrado (si no hay solución: "fuera del rango calibrado"), con

u^2(x_0) = \frac{\mathbf g^\top\mathbf C\,\mathbf g + s^2/m}{f'(x_0)^2},\qquad x_0 \pm t_{n-k}\,u(x_0)

( $\mathbf g = \partial f/\partial\boldsymbol\theta$ en $x_0$ , $m$ réplicas promediadas). Para la recta coincide exactamente con el $s_{x_0}$ clásico (validado con el ejemplo de Miller y Miller). Los campos LLOQ/ULOQ marcan los resultados fuera del intervalo de cuantificación.

Referencias: Findlay y Dillard, AAPS J. 9 (2007) E260; Miller y Miller, Statistics and Chemometrics for Analytical Chemistry, §5.6; Wald–Wolfowitz (prueba de rachas).

Cinética de Michaelis–Menten e inhibición

Código: matlens/core/bioassays.py::fit_michaelis_menten, fit_inhibition

v = \frac{V_\text{max}[S]}{K_m(1 + [I]/K_i) + [S](1 + [I]/(\alpha K_i))}

Competitiva ( $\alpha\to\infty$ ), acompetitiva (solo el término en $[S]$ ), no competitiva pura ( $\alpha = 1$ : $v = V_\text{max}[S]/\{(K_m+[S])(1+[I]/K_i)\}$ ) y mixta ( $\alpha$ libre). Ajuste global de todas las concentraciones de inhibidor (la concentración de inhibidor se lee de una columna de la hoja); $K_i$ con intervalo $t$ ; con "comparar" se ajustan los cuatro modelos y se elige el de menor AICc, mostrando los pesos de Akaike. Sin dos o más niveles de $[I]$ , $K_i$ no se calcula y se explica. Las representaciones Lineweaver–Burk, Eadie–Hofstee y Hanes–Woolf se dibujan con el ajuste no lineal (son solo diagnósticas).

Validación: el ajuste de Michaelis–Menten reproduce el resultado de Bates y Watts para los datos de puromicina (Treloar 1974).

Referencias: Copeland, Enzymes, 2.ª ed. (Wiley-VCH 2000), cap. 8; Bates y Watts, Nonlinear Regression Analysis and Its Applications (1988).

Curvas de crecimiento (Gompertz, logístico, Baranyi)

Código: matlens/core/bioassays.py::fit_growth

Con $y = \ln(N/N_0)$ (Zwietering):

y_\text{G} = A\exp\!\left\{-\exp\!\left[\frac{\mu_m e}{A}(\lambda - t) + 1\right]\right\},\qquad y_\text{L} = \frac{A}{1 + \exp\!\left[\frac{4\mu_m}{A}(\lambda - t) + 2\right]}

Con $y = \ln N$ (Baranyi y Roberts, $m=1$ , $h_0 = \mu_m\lambda$ ):

y = y_0 + \mu_m A(t) - \ln\!\left[1 + \frac{e^{\mu_m A(t)} - 1}{e^{y_\text{max}-y_0}}\right],\quad A(t) = t + \frac{1}{\mu_m}\ln\!\left(e^{-\mu_m t} + e^{-h_0} - e^{-\mu_m t - h_0}\right)

$\mu_m$ es la tasa específica máxima (pendiente máxima de $\ln N$ ), $\lambda$ la latencia (corte de la tangente en el punto de inflexión), $t_d = \ln 2/\mu_m$ con incertidumbre por el método delta. El diálogo toma $\ln$ de la señal (OD) salvo que se indique que ya está en escala logarítmica.

Referencias: Zwietering et al., Appl. Environ. Microbiol. 56 (1990) 1875; Baranyi y Roberts, Int. J. Food Microbiol. 23 (1994) 277.

Unión de ligandos (Kd, Bmax)

Código: matlens/core/bioassays.py::fit_binding

B_\text{total} = \frac{B_\text{max}[L]}{K_d + [L]} + NS\,[L]

Ajuste no lineal de un sitio; el término inespecífico $NS[L]$ es opcional. $K_d$ y $B_\text{max}$ con intervalos $t$ .

Referencia: Motulsky y Christopoulos (2004), cap. 42.

qPCR: Cq, eficiencia y expresión relativa

Código: matlens/core/bioassays.py::ct_threshold, ct_sdm, qpcr_efficiency, ddct, pfaffl

Umbral fijo: $C_t = c_{i-1} + (F_\text{th} - F_{i-1})(c_i - c_{i-1})/(F_i - F_{i-1})$ en el primer cruce.

Segunda derivada máxima (SDM): se ajusta $F = F_b + F_\text{max}/(1 + e^{-(c-c_0)/b})$ a toda la curva; el máximo de $F''$ está en

C_q^\text{SDM} = c_0 - b\,\ln(2+\sqrt3)

Eficiencia desde una curva estándar $C_q = a + m\log_{10}N_0$ (si no se da curva estándar no se muestra eficiencia y se explica por qué):

E = 10^{-1/m} - 1,\qquad u(E) = \frac{\ln 10\cdot 10^{-1/m}}{m^2}\,u(m)

$2^{-\Delta\Delta C_t}$ (supone $E = 100\,\%$ en ambos genes): $\Delta C_t = \bar C_{t,\text{diana}} - \bar C_{t,\text{ref}}$ con $s = \sqrt{s_d^2 + s_r^2}$ ; $\Delta\Delta C_t = \Delta C_{t,\text{muestra}} - \Delta C_{t,\text{calibrador}}$ conserva la $s$ de la muestra; intervalo $2^{-(\Delta\Delta C_t \pm s)}$ .

Pfaffl: $R = (1+E_d)^{\Delta C_{q,d}}/(1+E_r)^{\Delta C_{q,r}}$ con $\Delta C_q = C_q(\text{control}) - C_q(\text{muestra})$ , y

u^2(\ln R) = (\ln(1+E_d))^2u^2(\Delta C_{q,d}) + \Big(\frac{\Delta C_{q,d}}{1+E_d}\Big)^2u^2(E_d) + (\text{ídem referencia})

Referencias: Bustin et al. (MIQE), Clin. Chem. 55 (2009) 611; Livak y Schmittgen, Methods 25 (2001) 402; Pfaffl, Nucleic Acids Res. 29 (2001) e45; Spiess, Feig y Ritz, BMC Bioinformatics 9 (2008) 221; JCGM 100:2008 (GUM).

Algoritmo de ajuste en detalle

fit_nls, común a todos los modelos:

Los parámetros fijados se excluyen; se exigen más puntos que parámetros libres.
Los valores iniciales que tocan un límite finito se desplazan $10^{-12}(1+|\ell|)$ hacia dentro.
Sin límites se prueban Levenberg–Marquardt (lm) y región de confianza (trf) y se queda la solución de menor coste; con límites, solo trf. Tolerancias xtol = ftol = gtol $= 10^{-15}$ , hasta 20 000 evaluaciones, x_scale = "jac". Los residuos no finitos se sustituyen por $10^{150}$ .
En el óptimo el jacobiano se recalcula por diferencias centrales (paso $10^{-6}\max(|\theta|, 10^{-8})$ ) y $\mathbf C = s^2\,(\mathbf J^\top\mathbf J)^{+}$ (pseudoinversa).

Incertidumbre: qué incluye y qué no

Los errores estándar son asintóticos (linealización en el óptimo) con $s^2$ estimada de los residuos: suponen errores independientes, homocedásticos y de media nula. Los intervalos usan $t_{n-k}$ al nivel pedido (95 % por defecto), no $k = 2$ . No incluyen el sesgo del modelo (un 4PL sobre una curva asimétrica), la heterocedasticidad típica de las señales de placa (no hay ponderación) ni la incertidumbre de las concentraciones nominales. Las funciones derivadas (EC50 del 5PL, $t_d$ , $C_q$ SDM, $x_0$ ) usan el método delta, que puede quedarse corto con pocos grados de libertad. matlens/core/bioassays.py no emite indicadores warning/invalid; los avisos (prueba de rachas, fuera del rango calibrado, controles insuficientes) se muestran en los diálogos.

Ejemplo resuelto: dosis-respuesta 4PL

Curva sintética 4PL: $B = 5$ , $T = 100$ , $C_{50} = 0.5$ µM, $h = -1.2$ ; 11 concentraciones de $10^{-3}$ a $10^{2}$ µM (paso de media década) por triplicado, ruido gaussiano $\sigma = 3$ (semilla fija). fit_dose_response(..., '4pl'):

Parámetro	Valor	$u$ (1σ)	IC 95 %
bottom	5.18	1.24	[2.636, 7.724]
top	99.01	1.09	[96.78, 101.2]
logec50	-0.2995	0.0272	[-0.3551, -0.2438]
hill	-1.272	0.0893	[-1.455, -1.09]

EC50 = 0.5018 µM, IC 95 % [0.4415, 0.5704] µM (asimétrico); $R^2$ = 0.9937; $s$ = 3.45 con 29 grados de libertad.

Ejemplo resuelto: calidad de placa

Placa sintética: 8 controles positivos ~ N(1000, 40) y 8 negativos ~ N(150, 20) (semilla fija). plate_qc:

Magnitud	Valor
Z′	0.7234
SSMD	14.31
S/B	6.685
%CV (+)	5.43
%CV (−)	16.1
atípicos ( $\lvert M\rvert > 3.5$ )	0

Ejemplo resuelto: qPCR

qPCR: $C_t$ diana (muestra) 24.1, 24.3, 24.0; referencia (muestra) 18.0, 18.2, 17.9; diana (calibrador) 26.2, 26.0, 26.3; referencia (calibrador) 18.1, 18.0, 18.2. ddct:

Magnitud	Valor
$\Delta C_t$ muestra	6.1
$s$	0.216
$\Delta\Delta C_t$	-1.967
$2^{-\Delta\Delta C_t}$	3.909
intervalo $2^{-(\Delta\Delta C_t\pm s)}$	[3.365, 4.54]

Pfaffl con $E_d = 0.95 \pm 0.02$ , $\Delta C_{q,d} = 2.1 \pm 0.15$ , $E_r = 0.98$ , $\Delta C_{q,r} = 0$ : pfaffl da $R$ = 4.065, $u(\ln R)$ = 0.102, $u(R)$ = 0.417.