Commit c950aa82 by Turnhout, M.C. van

### start working it out on paper

parent 2ff6f85b
 function [amounts, imrgb] = clrdecon(dyes, imRGB) function [amounts, atoRGB, imrgb] = clrdecon(dyes, imRGB) % dye OD matrix: amount*nOD p = -log(dyes + 1); % add one to avoid taking log(0) (which would be bad) p = log(dyes + 1); % add one to avoid taking log(0) (which would be bad) % normalise: matrix M contains p-hat' % normalise: matrix M contains p-hat' or nOD' M = p; for r = 1:size(p, 2) % loop over rows (input colours) M(r, :) = M(r, :)/norm(p(r, :)); % divide row by its length end M = M % image OD imOD = -log(imRGB + 1); imOD = log(imRGB + 1); % dye contributions: 3 by 1 column % amounts = inv(M)*imOD amounts = M\imOD; % is faster and more accurate in Matlab amounts = inv(M')*imOD' amounts = M'\imOD'; % is faster and more accurate in Matlab % convert back to intensities atoRGB = zeros(3); for c = 1:3 % [amounts(1)*M(c, 1) amounts(3)*M(c, 3) amounts(3)*M(c, 3)] atoRGB(c, :) = exp(-amounts'.*M(c, :)) - 1; atoRGB(c, :) = exp(amounts'.*M(c, :)) - 1; % subtract one to scale back to 0-255 end imrgb = sum(atoRGB, 1); ... ...
 % !TeX root = colourdecon.tex %\begin{savequote} %The computable'' numbers may be described briefly as the real numbers whose expressions as a decimal are calculable by finite means.\qauthor{Alan Turing \cite{Turing1937}} %\end{savequote} \chapter{Linear algebra behind colour deconvolution}\label{algebra} \section{Introduction} Ruifrok and Johnston propose a method to do colour deconvolution for subtractive colour mixing in \cite{Ruifrok2001}. They use the law of Bouguer-Lambert-Beer \cite{Beer1852,Perrin1948} for a linear relationship between the amount dye' and the amount of absorbed light', and neglect such effects as reflection and scattering. \section{Transmittance and absorbance} The amount (fraction) of light that is transmitted through a certain amount of a certain dye $T$ is usually written as T(\lambda) = \frac{I(\lambda)}{I_0(\lambda)} = 10^{-\varepsilon(\lambda)cd} with $\lambda$\,[nm] the wavelength of the light, $I_0$\,[-] the intensity of the incoming light, $I$\,[-] the intensity of the light after the sample, $\varepsilon$\,[mol\ap{-1}] the molar extinction coefficient, $c$\,[mol/m] the dye concentration, and $d$\,[m] the path length through the sample. The corresponding absorbance' is then: A(\lambda) = \log\left(\frac{I(\lambda)}{I_0(\lambda)}\right) = -\log T(\lambda) = \varepsilon(\lambda)cd ~\\ \noindent We now choose to write this as: T(\lambda) = \mathrm{e}^{-k(\lambda)a} \label{Tlambda} With $k(\lambda)$\,[g\ap{-1}] the mass extinction coefficient, and $a$\,[g] the amount of dye in the sample (in grams), and with \begin{align} a & = cdM\\ k(\lambda) & =\frac{ \ln(10) \varepsilon(\lambda)}{M} \end{align} with $M$\,[g/mol] the molar mass of the dye, and $\ln(10)$ from the conversion of the base of the exponent. Absorbance is now linear in the extinction coefficient $k$ and the amount of dye $a$: A(\lambda) = -\ln T(\lambda) = k(\lambda)a \label{Alambda} In this formulation optical density' is an alternative term for $A(\lambda)$, and absorption coefficient' and attenuation coefficient' are alternative names for $k(\lambda)$. The concept that is exploited by Ruifrok and Johnson \cite{Ruifrok2001} is that: since optical density $A(\lambda)$ is linear in absorption $k(\lambda)$ and dye amount $a$ we can simple add contributions of $n$ different dyes to find their combined (total') absorbance or optical density. For $N$ dyes: A(\lambda) = \sum_{n=1}^N a_n k_n(\lambda) \section{Linear algebra with discrete wavelengths} For a discrete number of $N$ wavelengths, we can now write for the absorbance of a single dye in column format as: \begin{align} \col{A} & = \begin{bmatrix} A(\lambda_1)\\ A(\lambda_2)\\ \vdots \\ A(\lambda_N)\end{bmatrix} = \begin{bmatrix} A_\mathrm{1}\\ A_\mathrm{2}\\ \vdots \\ A_\mathrm{N}\end{bmatrix} \\ \col{k} & = \begin{bmatrix} k(\lambda_1)\\ k(\lambda_2)\\ \vdots \\ k(\lambda_N)\end{bmatrix} = \begin{bmatrix} k_\mathrm{1}\\ k_\mathrm{2}\\ \vdots \\ k_\mathrm{N}\end{bmatrix} \\ \col{A} = & a\col{k}\\ \end{align} However, this column for $\col{k}$ has a (definite) length in this discrete formulation, and this length affects the (strength of) absorption of a certain colour' as much as the factor $a$. For our purposes, we would like to separate dye amount' and 'dye colour effect' completely, i.e.\ only $a$ should affect the strength of absorption. So we divide $\col{k}$ by its length $l_k$ (normalize') to get $\col{\hat{k}}$ and add the original effect of this length to $a$ to get $\hat{a}$: \begin{align} l_k & = \sqrt{k_1^2 + k_2^2 + \dots +k_N^3} = \sqrt{\sum_{n=1}^N k_n^2}\\ \col{\hat{k}} & = \frac{\col{k}}{l_k} = \begin{bmatrix}\frac{k_1}{l_k}\\[.5em] \frac{k_2}{l_k}\\[.5em] \vdots\\[.5em] \frac{k_N}{l_k}\end{bmatrix} \\ \hat{a} & = al_k\\ \col{A} & = \hat{a}\hat{\col{k}} \end{align} Thus, the column $\hat{\col{k}}$ contains the \textsl{relative} absorption coefficients for our (three) discrete wavelengths such that the combined effect of the (three) absorption coefficients adds up to 1. The total absorption $\col{A}$ for our $N$ discrete wavelengths for $D$ different dyes can now be written as \col{A} = \sum_{d=1}^D\hat{a}_d \col{\hat{k}}_d \label{sumdyedis} Which simply says that the total amount of light absorbed at a certain wavelength is the sum of the amount absorbed of the first dye at that wavelength plus the amount absorbed of the second dye at that wavelength plus the amount absorbed of the third dye at that wavelength plus \dots\\ Equation \ref{sumdyedis} can be written in matrix form as: \begin{align} \mat{K} & = \begin{bmatrix} \col{\hat{k}}_1 & \col{\hat{k}}_ 2 & \dots & \col{\hat{k}}_D\end{bmatrix},\; \col{a} = \begin{bmatrix} \hat{a}_1 \\ \hat{a}_ 2\\ \vdots \\ \hat{a}_D\end{bmatrix}\\ \mat{K}\col{\hat{a}} & = \col{A} \label{sumdyesmat} \end{align} Matrix $\mat{K}$ is composed of the $D$ columns with (normalized) absorption coefficients $\col{\hat{k}}$ put together side-by-side. Since each column $\col{\hat{k}}$ contains $N$ coefficients for $N$ wavelengths, the size of $\mat{K}$ is $N\times D$. The column $\hat{a}$ contains the $D$ dye amounts $\hat{a}_1$ through $\hat{a}_D$, and column $\col{A}$ contains the total (summed) absorption of the $D$ dyes for each of the $N$ wavelengths. \noindent We can always calculate total absorbance $\col{A}$ for known $\mat{K}$ and $\hat{a}$ with equation \ref{sumdyesmat}, but to find the amounts $\hat{a}$ for given $\col{A}$ and $\mat{K}$, we have to solve: \begin{align} \mat{K}\col{\hat{a}} & = \col{A} \tag{\ref{sumdyesmat}}\\ \mat{K}^{-1}\mat{K}\col{\hat{a}} & = \mat{K}^{-1}\col{A} \\ \col{\hat{a}} & = \mat{K}^{-1}\col{A} \label{sumdyesinv} \end{align} and the inverse $\mat{K}^{-1}$ only exists for square matrices $\mat{K}$. Not to worry, we can get a square $\mat{K}$ by simply evaluating as much dyes $D$ as wavelengths $N$: $D = N$. \section{Colour deconvolution in RGB} \subsection{Deconvolution in absorption-space'} When we use three wavelengths and label them, say, $R$, $G$, and $B$ and (thus) use three dyes labelled 1, 2 and 3, equation \ref{sumdyesmat} comes out as: \begin{bmatrix} \hat{k}_{R_1} & \hat{k}_{R_2} & k_{R_3} \\ \hat{k}_{G_1} & \hat{k}_{G_2} & \hat{k}_{G_3} \\ \hat{k}_{B_1} & \hat{k}_{B_2} & \hat{k}_{B_3}\end{bmatrix}\cdot\begin{bmatrix}\hat{a}_1\\\hat{a}_2\\\hat{a}_3 \end{bmatrix} = \begin{bmatrix}\hat{a}_1 \hat{k}_{R_1} + \hat{a}_2\hat{k}_{R_2} +\hat{a}_3 \hat{k}_{R_3} \\\hat{a}_1 \hat{k}_{G_1} + \hat{a}_2\hat{k}_{G_2} +\hat{a}_3 \hat{k}_{G_3} \\ \hat{a}_1 \hat{k}_{B_1} + \hat{a}_2\hat{k}_{B_2} +\hat{a}_3 \hat{k}_{B_3} \end{bmatrix}= \begin{bmatrix}A_R\\A_G\\A_B \end{bmatrix} and equation \ref{sumdyesinv} comes out as: \begin{bmatrix}\hat{a}_1\\\hat{a}_2\\\hat{a}_3 \end{bmatrix} = \begin{bmatrix} \hat{k}_{R_1} & \hat{k}_{R_2} & k_{R_3} \\ \hat{k}_{G_1} & \hat{k}_{G_2} & \hat{k}_{G_3} \\ \hat{k}_{B_1} & \hat{k}_{B_2} & \hat{k}_{B_3}\end{bmatrix}^{-1}\cdot \begin{bmatrix}A_R\\A_G\\A_B \end{bmatrix} So when we know matrix $\mat{K}$ (and can calculate its inverse), and know the total amount of absorption for the three wavelengths $\col{A}$, we can calculate the contributions (amounts) of the three dyes $\hat{a}$. Q.E.D.! \subsection{Except that our images are in transmission-space'} Except that your RGB images are in transmission space': the RGB scheme uses additive mixing. Worse, it is fundamentally impossible to measure absorption: you always measure intensity, brightness, transmission', not the amount of darkness'.\\ \noindent We can use equations \ref{Tlambda} and \ref{Alambda} to go from transmission-space' to absorption-space' and back. When the maximum possible pixel intensity is $I_\mathrm{max}$ ($2^n-1$ for a $n$-bits camera) and the actual/recorded pixel intensity (RGB) column is $\col{I}_p$, then: \col{T} = \frac{\col{I}_p}{I_\mathrm{max}} = \mathrm{e}^{-\col{\hat{k}}\hat{a}} which would make pixel intensity from absorption \col{I}_p = I_\mathrm{max} \col{T} = I_\mathrm{max} \mathrm{e}^{-\col{\hat{k}}\hat{a}} and absorption from pixel intensity: \col{A}= -\ln \left(\frac{\col{I}_p}{I_\mathrm{max}} \right) = \ln I_\mathrm{max} -\ln \col{I}_p =\col{\hat{k}}\hat{a} ~\\ Except that pixel values can be zero, and the (any) log of zero is $-\infty$ which makes numerical computations hard. So we add 1 to the pixel intensities to avoid that we have to store $\pm\infty$ in the computer: \col{A}= -\ln \left(\frac{\col{I}_p + 1}{I_\mathrm{max}} \right) = \ln I_\mathrm{max} - \ln \left(\col{I}_p+ 1\right) \approx \col{\hat{k}}\hat{a} But \ No newline at end of file
 % Encoding: UTF-8 @Article{Beer1852, author = {Beer, August}, title = {{B}estimmung der {A}bsorption des rothen {L}ichts in farbigen {F}l{\"u}ssigkeiten}, journal = {Annalen der Physik}, year = {1852}, volume = {162}, number = {5}, pages = {78--88}, doi = {10.1002/andp.18521620505}, owner = {mark}, timestamp = {2010.06.07}, } @Article{Ruifrok2001, author = {Ruifrok, Arnout C. and Johnston, Dennis A.}, title = {Quantification of histochemical staining by color deconvolution}, journal = {Analytical and quantitative cytology and histology}, year = {2001}, volume = {23}, number = {4}, pages = {291--299}, abstract = {Objective: To develop a flexible method of separation and quantification of immunohistochemical staining by means of color image analysis. Study Design: An algorithm was developed to deconvolve the color information acquired with RGB cameras, to calculate the contribution of each of the applied stains, based on the stain-specific RGB absorption. The algorithm was tested using different combinations of DAB, hematoxylin and eosin at different staining levels. Results: Quantification of the different stains was not significantly influenced by the combination of multiple stains in a single sample. The color deconvolution algorithm resulted in comparable quantification independent of the stain combinations, as long as the histochemical procedures did not influence the amount of stain in the sample due to bleaching because of stain solubility, and saturation of staining was prevented. Conclusion: The presented image analysis algorithm provides a robust and flexible method for objective immunohistochemical analysis of samples stained with up to three different stains, using a laboratory microscope and standard RGB camera setup, and the public domain program NIH image.}, keywords = {image processing, computer-assisted; color separation, color deconvolution, immunohistochemistry}, owner = {tue}, pmid = {11531144}, publisher = {St. Louis, MO: Science Printers and Publishers, 1985-}, timestamp = {2020.09.18}, } @Article{Perrin1948, author = {Perrin, Fred H.}, title = {{W}hose {A}bsorption {L}aw?}, journal = {Journal of the Optical Society of America}, year = {1948}, volume = {38}, number = {1}, pages = {72--74}, doi = {10.1364/JOSA.38.000072}, owner = {mark}, timestamp = {2010.06.07}, } @Misc{Oransky2017, author = {Oransky, Ivan and {Marcus (\href{https://retractionwatch.com}{Retraction Watch})}, Adam}, title = {Welcome to the {J}ournal of {A}lternative {F}acts. {T}hey're the greatest! {A}nd winning!}, howpublished = {\newline \href{https://retractionwatch.com/2017/01/31/welcome-journal-alternative-facts-theyre-greatest-winning/}{https://retractionwatch.com/2017/01/31/welcome-journal-alternative-facts-theyre-greatest-winning/}}, month = {January}, year = {2017}, owner = {tue}, timestamp = {2020-08-16}, } @Comment{jabref-meta: databaseType:bibtex;}