\section{Introduction}
\vspace{-0.2cm}
Breast cancer is the most commonly diagnosed cancer for women, which is estimated to account for 30\% of new cancer diagnoses and 15\% of cancer deaths in the United States \cite{RN4}. Early and precise diagnosis of breast cancer is crucial to improve the survival rate of patients \cite{RN2}. Microscopic examination of stained tissue sections is among the most accurate methods of diagnosing and classifying cancer.
The cancer effects can be observed in WSIs in the cellular and tissue levels. 
\begin{figure}
\centering
\includegraphics[width=0.8\textwidth]{intro1.PNG}
\setlength{\abovecaptionskip}{-0.2cm} 
\caption{An annotated whole slide image: The patches framed in red contour are invasive cancer. Those framed in blue contour are  
In-Situ cancer. Those framed in green contour are Benign. 
The rest part of the slide is normal.} \label{fig1}
\vspace{-0.6cm}
\end{figure}

Fig.\ref{fig1} shows an example of the cancer effects, which are classified into four categories, i.e., invasive cancer, In-Situ cancer, Benign, and normal. 
Recently, many computer-aided systems utilize deep learning models to improve the classification consistency and accuracy \cite{RA9_MT,RA10_MT,RA16_MT,RA9_MT}. However, robust deep learning models require large annotated datasets, which are costly to produce especially for medical images. Recent studies integrate active learning with deep learning, which utilize unannotated data to improve the performance of deep learning model \cite{RA1_AL,RA2_AL,RA11_AL,RA13_AL}. 
Yang et al. applied active learning method on fully convolutional network (FCN) to select the most representative and uncertain areas for annotation \cite{RA6_FCN}. One of the drawbacks is that the FCN cannot be applied to high resolution images, such as WSIs. And it is difficult to acquire iterative annotation on the high-resolution WSIs. Generative Adversarial Network (GAN) is widely applied to generate realistic images, which overcomes the limitations of small training datasets. Mahapatra et al. applies conditional generative adversarial networks (cGANs) to generate informative and realistic chest X-ray images, which enlarge the training datasets \cite{RA3_GAN}. 
However, GAN generates artifacts on the original images. The artifacts introduced to the images, although subtle, could substantially alter the features of cells and tissues, which mislead the model and affect the convergence of parameters. 
Existing methods fail to provide an efficient solution for automatically grading breast cancer on limited annotated WSIs.
In this paper, we present a new semi-supervised framework based on expectation maximization (EM). We leverage unannotated WSIs to adjust the deep learning model on a limited annotated dataset, reducing the reliance on expensive pixel-wise annotations.
The main contributions of the paper are: 
\begin{itemize}
\item To the best of our knowledge, we first apply active learning method in breast cancer grading, and propose a semi-supervised approach based on EM to effectively reduce annotated dataset for multi-classes pixel-wise breast cancer grading on WSIs.
\item We propose a sample selection method based on collaborative filtering, which selects the credible and representative unannotated datasets for enlarging the training dataset.
\item Using the proposed semi-supervise framework, significantly enhanced performance on pixel-wise prediction of WSIs is achieved with only 30\% of annotated dataset. 
\end{itemize}

\begin{figure}
\includegraphics[width=\textwidth]{method2.PNG}
\setlength{\abovecaptionskip}{-0.4cm} 
\caption{Overall flow of the EM model.} 
\label{Figure.2}
\vspace{-0.6cm}
\end{figure}


\vspace{-0.5cm}
\section{Method}
\vspace{-0.4cm}
\subsection{Semi-supervised Learning Framework Based on EM model}
\vspace{-0.2cm}
An overview of our semi-supervised learning framework is shown in Fig.\ref{Figure.2}  
In the semi-supervised learning framework, only part of whole slide images is annotated, which is defined as set $D$. The label of some slides is unknown, which is defined as set $U$. 
Let $y_i$ denote the label for patch $x_i \in D$. 
Let hidden variable $z_j$ denote the label for patch $x_j \in U$. 
We initialize the CNN model on $D$ and update the model parameter to $\theta^0$. We apply initial CNN model to produce the probability map $P(z_j|x_j)$ of $x_j \in U$. 
The EM algorithm alternates between the E-step for estimating the hidden labels $z_j$ and the M-step for computing optimal model parameters with maximized $P(X|\theta, Z)$.
The probability map $P(z_j|x_j)$ is projected to a scaled value between 0 and 1, which is used to generate the consistent heatmap (see Fig.\ref{Figure.2}(C)).
The fixed vector $\beta^*$ = $(\beta^1, \beta^2, \beta^3)$ is applied on the heatmap to generate the classmap as shown in Fig.\ref{Figure.2}(D).

Next, The most representative and credible patches based on collaborative filtering are selected to train the CNN model in the next iteration.\\
\textbf{Initialization}
Assume the patches are independently and identically distributed (i.i.d.). The initial parameter $\theta^0$ is obtained from the CNN model, which is trained on annotated dataset $D$. Here, $\theta^0$ is computed as: 
\begin{equation}
\begin{aligned}
\theta^0 & \leftarrow  \mathop{\arg\max}_{\theta} \prod_{x_i \in D}P(x_i,y_i| \theta)
&= \mathop{\arg\max}_{\theta} \prod_{x_i \in D}P(y_i| x_i;\theta)P(x_i|\theta)
\end{aligned}
\label{eqn:theta}
\end{equation}
\textbf{E-step:} Based on the current parameters $\theta^t$ at EM iteration $t$, we calculate the probability maps $P(z_j | x_j, \theta^t )$ of unannotated patches, and then re-scale the probability maps to $P_{norm}(z_j | x_j, \theta^t) \in [0,1]$. We generate the class label $c(x_j)$ based on $\beta^*$ and then obtain the classmap.

The ground truth of the unannotated patches is then extracted as $c^j$, and the effective dataset $E^t$ is selected according to the method described in Section 2.2.\\
\textbf{M-step:} 
The CNN model is retrained on the effective dataset $E^t$ produced in the E-step. The model parameter $\theta$ is updated to maximize the likelihood defined in Equation (\ref{eq:M}). 
\begin{equation}
\begin{aligned}
\theta^{t+1}& = \mathop{\arg\max}_{\theta} Q(\theta, \theta^t)\\
K(\theta, \theta^t) &=  \prod_{x_i \in D}P(x_i,y_i| \theta) \times \prod_{x_j \in E^t}P(x_j,z_j| \theta)\\
&=  \prod_{x_i \in D}P(y_i|x_i; \theta)P(x_i|\theta) \times \prod_{x_j \in E^t}P(z_j| x_j; \theta)P(x_j|\theta)
\end{aligned}
\label{eq:M}
\end{equation}
Assume that $x_j | \theta$ follows an uniform distribution, we formulate the objective function $Q(\theta, \theta^t)$ as:
\begin{equation}
\begin{aligned}
    Q(\theta, \theta^t) &=\log K(\theta, \theta^t)\\
    & \propto \sum_{x_i \in D}\log P(y_i|xi;\theta) + \sum_{x_j \in E^t}\log P(z_j|x_j;\theta)
\end{aligned} 
\label{eq:Q}
\end{equation}


\subsection{Patch Selection}
\vspace{-0.2cm}
\textbf{Hard Example Mining}
Hard example mining is used in the initialization step to fully exploit the annotated dataset, especially those with wrong classification results.
An effective coefficient $\alpha$ is defined as in Equation (\ref{eq:A}). The higher the value of $\alpha$ is, the harder and more valuable the corresponding patch is for model training.
\begin{equation}
\begin{aligned}
& c_k = \mathop{\arg\max}_{j} p(y_i^{c_j})
& \alpha =  \|c_k - c_i\| \times P(y_i^{c_k})
\end{aligned}
\label{eq:A}
\end{equation}
Here, $c_i$ denotes the class label of the patch $x_i$, and $P(y_i^{c_j})$ denotes the probability map.

For the initialization step, we first train our model on 50\% of the annotated data. Then, we apply this model on the rest of the data and calculate the effective coefficients. Patches with effective coefficient in the first quintile (top 20\%) are selected to retrain the model.\\
\textbf{Collaborative Filtering}
\begin{algorithm}[!h]
  \caption{Patch selection method.}
  \begin{algorithmic}[1]
    \renewcommand{\algorithmicrequire}{\textbf{Input:}}
    \REQUIRE ~\\
      $\mathcal{U}=\{x_i\},i\in[1,n]$~(Unannotated dataset)\\
      $\mathcal{D}=\{(y_j,L_j)\},j\in[1,m]$~(Annotated dataset)\\
      $\mathcal{M}_t$~(CNN model in iteration $t$)\\
      $\sigma$~(Similarity threshold for patch selection)
    \renewcommand{\algorithmicrequire}{\textbf{Output:}}
    \REQUIRE ~\\
      $\mathcal{E}_t$~(Set of unannotated patches in iteration $t$)
    \renewcommand{\algorithmicrequire}{\textbf{Functions:}}
    \REQUIRE ~\\
      $feature\gets F(\mathcal{M}, x)$~\{Output $512\times1\times1$ feature of $\mathcal{M}$ given patch $x$\}\\
      $prediction\gets P(\mathcal{M}, x)$~\{Prediction result of $\mathcal{M}$ for patch $x$\}\\
      $s\gets\mathrm{sim}(\mathbf{x},\mathbf{y})~\{s=\frac{\mathbf{x}\cdot\mathbf{y}}{\left\|\mathbf{x}\right\|\times\left\|\mathbf{y}\right\|}\}$\\
      $label\gets\mathrm{argmax_{index}}(num)$~\{Output label with largest number\}
    \renewcommand{\algorithmicrequire}{\textbf{Initialize:}}
    \REQUIRE ~\\
      $\mathcal{E}^t\gets\emptyset$
    \FOR {each $x_i\in\mathcal{U}$}
      \STATE $\alpha_i\gets F(\mathcal{M}^t,x_i)$
      \STATE $pred_i\gets P(\mathcal{M}^t,x_i)$
     
      \STATE Set $num$ to vector [0,0,0,0]
      \FOR {each $(y_j,L_j)\in\mathcal{D}$}
        \STATE $\beta_j\gets F(\mathcal{M}^t,y_j)$
        \IF {$\mathrm{sim}(\alpha_i,\beta_j)>\sigma$}
          \STATE $num[L_j]\gets num[L_j]+1$
        \ENDIF
      \ENDFOR
      \STATE $label_i\gets argmax(num)$
      \IF {$pred_i=label_i$}
        \STATE $\mathcal{E}^t\gets \mathcal{E}^t\cup x_i$
      \ENDIF
    \ENDFOR
	\end{algorithmic}
\end{algorithm}
In the E-step, patches are selected using Algorithm 1. We first apply CNN to extract the features of all patches, and then calculate similarity $sim(x_i,y_j)$ between each unannotated patch $x_i$ and annotated patch $y_j$. For each unannotated patch $x_i$, we compute the set of annotated patches as $\{y_j|sim_{y_j \in D}(x_i,y_j) >t\}$. Then we apply the majority voting method on the above computed patches to determine the label of the unannotated patch. If the assigned label is consistent with the predicted one by the model, we insert unannotated patch $x_i$ into $E^t$.

\vspace{-0.5cm}



\section{Experiment}
\vspace{-0.2cm}
\subsection{Data Preprocess}
\vspace{-0.2cm}
\textbf{Foreground Patch Extraction on WSI Image} The high resolution images of WSIs need to be converted into patches for use.
However, a large part of WSI is background, which produces uninformative patches in the datasets, and thus should be excluded. 
The widely used foreground extraction method Ostu fails to extract certain parts of the tissue from the slide. 
Therefore, we adapt the graph-based image segmentation method in \cite{RA12_FE} for foreground extraction. 
For a given slide, we construct an undirected graph $G(V,E)$. In $G$, each node $v_{i,j} \in V$ corresponds to a pixel. The edge set $E = \{(v_{i,j},v_{i+1,j}), (v_{i,j}, v_{i,j+1}))\}$ correspond to the connection between adjacent pixels. We set the edge weight to be $W(v_{i,j}, v_{i+1,j}) = \|v_{i,j} - v_{i+1,j}\|$.
Then we compute the minimum spanning tree $T$ using Kruskal's algorithm \cite{RA14_FE}, and delete the edges in $T$ whose weights are greater than a prespecified threshold (100 in the experiments). 
The deletion of these edges produces a forest, i.e., a set of sub-trees (e.g., $T_1, T_2, \cdots, T_n$). Now we compute the average RGB values for the sub-trees (e.g., $RGB(T_1), RGB(T_2), \cdots, RGB(T_n)$). Among the computed average RGB values, assume the maximum value is $u$. Then all the sub-trees with average RGB value greater than $u-45$ are set as background. Then the foreground mask is obtained as shown in Fig.\ref{Figure.4}.
According to the foreground mask, we crop the WSI into patches with 50\% overlap, where each patch consists of $1536 \times 1536$ pixels.
Patches with less than 40\% foreground pixels are considered to be background, which are not used for classification.\\
\begin{figure}[ht]
\centering
\includegraphics[width=\textwidth]{mask.PNG}
\caption{The foreground extraction of a WSI: (a) The original WSI, (b) Otsu method \cite{RA15_FE}, and (c) our method.} \label{Figure.4}
\vspace{-0.4cm}
\end{figure}
\textbf{Patch Label Extraction}
We assign the label of the patches according to the ground-truth contour of WSI. In most cases, the label of the patch is obtained according to the type of cancer (Benign, In-Situ or Invasive) with the largest area in the patch. However, there are two special cases as follows:
\begin{itemize}
\item If the cancer area in a patch is less than one-third of the whole patch area, this patch is labeled as normal.
\item If there are two or more types of cancer in a patch, and the corresponding tissue areas are both greater than one-third of the whole patch area, this patch is considered to be a noisy patch and discarded. In the experiments, the number of such patches is very few. 
\end{itemize}

In Kwok's work in \cite{RA9_MT}, the class value of a patch is the mean of the class values of all pixels in the patch. Our experiments show that the above method tends to generate wrong labels, which disrupt the learning process. For example, when half of a patch contains In-Situ areas and the rest is normal, Kwok's method labels this patch as benign even if there are no benign tissues at all.
In contrast, our EM-based method effectively avoids the drawbacks of Kwok's method.

\renewcommand{\arraystretch}{1.5}
\begin{table}[tp]
  \centering
  \fontsize{10}{8}\selectfont
  \begin{threeparttable}
  \caption{Demographic Prediction performance comparison by three evaluation metrics.}
  \label{tab:performance_comparison}
    \begin{tabular}{|c|c|c|c|c|}
    \toprule
    \multirow{1}{*}{\raisebox{-1.5ex}[0pt]{Metric}} &
    \multicolumn{2}{c|}{Patch-wise} & \multicolumn{2}{c|}{Pixel-wise} \cr
    \cmidrule(lr){2-5}
    & Precision & Accuracy & F1-Measure & Score Metric\cr
    \midrule
    Kwok et al. &0.5798&0.7784&0.6391&0.7715\cr
    ALM-10 \% &0.7205&0.8466&0.8078&0.7186\cr
    ALM-20 \% &0.7350&0.8684&0.8082&0.7447\cr
    ALM-30 \% &0.7477&0.9035&0.8303&0.7759\cr
    FSL &0.8239&0.8963&0.7698&0.7592\cr
    Our-10 \%&0.7218&0.8856&0.8054& 0.7675\cr
    Our-20 \% &0.7499&0.8852&0.8048& 0.7539\cr
    Our-30 \%&{\bf 0.7987}&{\bf 0.9197}&{\bf 0.8623}&{\bf 0.7858}\cr
    iteration2-Our-30\% &0.8293&0.9210&0.8751&0.8027\cr
    \bottomrule
    \end{tabular}
    \end{threeparttable}
\vspace{-0.6cm}
\end{table}

\vspace{-0.6cm} 
\begin{figure}
\includegraphics[width=\textwidth]{result.PNG}
\caption{The pixel-wise classification results of slides A02 and A08 (Green contour: Benign, Red contour: Invasive, Blue contour: In-Situ, Others: Normal), (a) and (e) give the results of Kwok methods, (b) and (f) give the results of the model trained on full annotated datasets, (c) and (g) give the results of our EM-based method trained on 30\% of annotated datasets combined with unannotated datasets, (d) and (h) give the results labeled by pathologist.} \label{fig.6}
\vspace{-0.6cm}
\end{figure}


\subsection{Patch Classification}
\vspace{-0.2cm}
Our neural network is a fine-tuned vgg19 network with batch normalization. 
Given the large patches of size $1536\times1536$, we resize them into $512\times512$, and then feed them into the network. 
For adapting to the fully connected layers, we add an average pooling layer, which converts the $512\times16\times16$ feature map into a $512\times1\times1$ vector.
The patch-wise experimental results are summarized in Tabel 1. 
we first apply active learning method (ALM) to continually finetune the classification model with informative and effective datasets instead of retrain the model with all datasets.
ALM-10\%, ALM-20\%, and ALM-30\% refer to different models trained on corresponding portions of annotated datasets. FSL is a model trained on all the annotated dataset. Our-10\%, Our-20\%, and Our-30\% are the proposed EM-based model trained on corresponding portions of annotated datasets. For example, Our-30\% denotes our proposed method using 30\% of the whole annotated dataset.

From the experiment, FSL obtains 0.76, 0.89 and 0.82 for F1 score, accuracy and precision, respectively. 
Kwok's method obtains similar results of 0.63, 0.77, and 0.57.
However, with 30\% of the whole annoatated dataset, ALM obtains 0.83, 0.90, and 0.74, which outperforms FSL. This can be explained by the exclusion of the uninformative data. 
In contrast, our proposed method obtains 0.86, 0.91, 0.79 for F1 score, accuracy and precision, respectively. 
Among the different methods, Our-30\% achieved the best results using only 30\% of annotated data combined with unannotated data. Moreover, our method significantly reduces the runtime for finetuning the model as in ALM.

\vspace{-0.4cm}
\subsection{Pixel-wise Classification on WSI}
\vspace{-0.2cm}
In the experiment, we first extract the probability map from the softmax layer of our model, and then normalize the probability map to a scaled value. The scaled value represents the intensity of cancer, where the scaled value of Normal cancer is close to 0, and the scaled value of Invasive cancer is close to 1. Given the scaled value of each patch, we generate a heatmap of WSI, where the intensity of every pixel is the scaled value of the patch that the pixel belongs to. Next, we map the heatmap to classmap by the fixed vector $\beta^* = (0.1, 0.5, 0.75)$, where a pixel is classified according to its intensity value as follows: (1) $[0,0.1]$ for Normal, (2) $(0.1,0.5]$ for Benign, (3) $(0.5,0.75]$ for In-Situ, and $(0.75,1]$ for Invasive.
The pixel-wise classification results are summarized in Table 1 in terms of score metrics defined on ICIAR. Our method with 30\% annotated dataset achieves the best performance with a score of 0.785, where the best score of Kwok, FSL, and ALM methods are 0.771,0.759 and 0.775, respectively.

The three methods Kwok, FSL and ALM perform relatively well in detecting large areas of cancer. However, for small areas of cancer, these methods usually fail. Fig.\ref{fig.6} shows an example of invasive tissues in A08 slide, which consists of many small cancer areas. Kwok's method tends to classify small invasive tissues to In-Situ tissues. On the other hand, FSL is unable to recognize lots of small invasive tissues in A02 and A08 slides shown in Fig.\ref{Figure.4}. In contrast, the proposed EM-based method is able to detect small areas of cancer, which are crucial for correct diagnosis.

\vspace{-0.5cm}
\section{Conclusion}
\vspace{-0.2cm}
In this paper, we have proposed an effective semi-supervised approach based on the EM model, which significantly reduces the reliance on the annotated dataset. Experiment results show that the proposed method achieves remarkable performance with only 30\% annotated datasets. 
Moreover, the proposed method effectively traces the small cancer areas, which is one of the key markers for cancer diagnosis. 
In the future, more parameters and metrics will be introduced in the system, such as max area of cancer, number of different types of cancer, degree of patient, etc. More prior knowledge will be introduced for generating adaptive parameters in the proposed EM framework.

\vspace{-0.5cm}
\bibliographystyle{splncs04}
