# SOSPCNN: Structurally Optimized Stochastic Pooling Convolutional Neural Network for Tetralogy of Fallot Recognition

Authors: Shui-Hua Wang; Kaihong Wu; Tianshu Chu; Steven L. Fernandes; Qinghua Zhou; Yu-Dong Zhang; Jian Sun
Journal: Wireless Communications and Mobile Computing
 (2021)

---

## Abstract

Aim. This study proposes a new artificial intelligence model based on cardiovascular computed tomography for more efficient and precise recognition of Tetralogy of Fallot (TOF). Methods. Our model is a structurally optimized stochastic pooling convolutional neural network (SOSPCNN), which combines stochastic pooling, structural optimization, and convolutional neural network. In addition, multiple-way data augmentation is used to overcome overfitting. Grad-CAM is employed to provide explainability to the proposed SOSPCNN model. Meanwhile, both desktop and web apps are developed based on this SOSPCNN model. Results. The results on ten runs of 10-fold crossvalidation show that our SOSPCNN model yields a sensitivity of 92.25±2.19, a specificity of 92.75±2.49, a precision of 92.79±2.29, an accuracy of 92.50±1.18, an F1 score of 92.48±1.17, an MCC of 85.06±2.38, an FMI of 92.50±1.17, and an AUC of 0.9587. Conclusion. The SOSPCNN method performed better than three state-of-the-art TOF recognition approaches.


---

## Body

## 1. Introduction

Tetralogy of Fallot (TOF) is a congenital defect that influences normal blood flow through the heart [1]. It is made up of 4 defects of the heart and its blood vessels [2]: (a) ventricular septal defect, (b) overriding aorta, (c) right ventricular outflow tract stenosis, and (d) right ventricular hypertrophy. Defects of TOF can cause oxygen in the blood that flows to the rest of the body to be reduced. Infants with TOF have a bluish-looking skin color [3] since their blood does not carry enough oxygen.Traditional diagnosis of TOF is after a baby is born, often after the infant had an episode of cyanosis during crying or feeding. The most common test is an echocardiogram [4], an ultrasound of the heart that can show problems with the heart structure and how well the heart is working with this defect. Recently, computed tomography (CT) has shown its success in the differential diagnosis of TOF [5], since it can provide detailed images of many types of cardiovascular issue; besides, computed tomography (CT) can be performed even if the subject has an implanted medical device, unlike magnetic resonance imaging (MRI) [6].Manual diagnosis on CT is lab-intensive, onerous, and needs expert skills. Besides, the manual results vary due to intraexpert and interexpert factors. Shan et al. (2021) [7] mention that “fully manual delineation that often takes hours” and the modern automatic diagnosis models based on artificial intelligence (AI) can only take seconds to minutes to get decisions, which now becomes a hot research field.For example, Ye et al. (2011) [8] present a morphological classification (MC) method. The authors extract morphological features by registering cardiac MRI scans to a template. Later, deep learning (DL) rises as a new type of artificial intelligence (AI) technique and has shown its powerfulness in many academic and industrial fields. Within the field of DL, convolutional neural network (CNN) is one standard DL algorithm that is particularly suitable for handling images. Giannakidis et al. (2016) [9] presented a multiscale three-dimensional CNN (3DCNN) for segmentation of the right ventricle. Tandon et al. (2021) [10] present a ventricular contouring CNN (VCCNN) algorithm.The difference between this study to previous studies is that we simplify the problem to a binary-coded classification problem [11]; that is, given an input cardiovascular CT image, the AI model should have the ability to give a binary output, i.e., predict whether the subject is TOF or healthy. This simplification makes the AI model focus on the prediction task itself and does not need to generate human-understandable outputs (such as segmentation, contouring, etc.) in the light of the expectation to make our AI model more accurate. Furthermore, we propose a new stochastic pooling CNN (SCCNN) that uses a new pooling technique—stochastic pooling to improve the prediction performance. All in all, our contributes are fourfold:
(a)
Stochastic pooling is employed to replace traditional max-pooling(b)
Structural optimization is carried out to fix the optimal structure(c)
Multiple-way DA is introduced to increase the diversity of training images(d)
Experiments by ten runs of 10-fold crossvalidation show that our method is better than three state-of-the-art approachesThe rest of this paper is structured as follows: Section2 describes the dataset. Section 3 contains the rationale of methodology, including the preprocessing, stochastic pooling, structural optimization, multiple-way data augmentation, the implementation, Grad-CAM, and evaluation measures. Section 4 presents the experimental results and discussions. Section 5 concludes this paper.


## 2. Dataset

This study is a retrospective research, of which ethical approval is exempted. The imaging protocol is described below: Philips Brilliance 256 row spiral CT machine, KV: 80, MAS: 138, Layer Thickness 0.8 mm, Lung Window (W: 1600 HU, L: -600 HU), Mediastinal Window (W: 750 HU, L: 90 HU), thin layer reconstruction according to the lesion display, layer thickness, and layer distance are both 0.8 mm mediastinal window images. Place the patient in a supine position, let the patient breathe deeply after holding in, and conventionally scan from the apex of the lung to the costal diaphragmatic angle. The resolutions of all images are 512 by 512 pixels. Data is available upon reasonable requests to corresponding authors.We selected ten children with Tetralogy of Fallot who were admitted to Nanjing Children’s Hospital from March 2017 to March 2020. We then used a systematic random sampling method to select ten normal children from healthy medical examiners within the same period of time. The Tetralogy of Fallot (TOF) observation group included three males and seven females, aged 4-22 months, with an average age of (8.90±5.47) months. Normal children in the control group included six males and four females, aged from 3 months to 24 months, with an average age of 10.4±8.14 months. Inclusion criteria for children with confirmed Tetralogy of Fallot are as follows:
(1)
CT suggests Tetralogy of Fallot(2)
Surgery confirmed that the anatomical deformity of the heart is Tetralogy of Fallot


## 3. Methodology

### 3.1. Preprocessing

Table1 lists the abbreviation list for the ease of reading. A five-step preprocessing was carried out on all the images to select the important slices, save storage, enhance contrast, remove unnecessary image regions, and reduce the image resolution.Table 1
Abbreviation list.
AbbreviationMeaningAPAverage poolingAUCArea under the curveBNNormalizationCNNConvolutional neural networkCTComputed tomographyDAData augmentationDPDDiscrete probability distributionFCLFully connected layerFMIFowlkes–Mallows indexGrad-CAMGradient-weighted class activation mappingGUIGraphical user interfaceHSHistogram stretchingL2Pl2-norm poolingMCCMatthews correlation coefficientMPMax-poolingMRIMagnetic resonance imagingMSDMean and standard deviationPMProbability mapReLURectified linear unitRLVRandom location vectorROCReceiver operating characteristicSCStrided convolutionSPStochastic poolingTOFTetralogy of FallotFirst, four slices were chosen by radiologists using a slice-level selection method. For TOF patients, the slices showing the largest size and number of lesions were selected. For healthy control subjects, any level of the image can be selected. Now, we have in total 40 TOF images and 40 HC images.Second, all the images are converted to grayscale images and stored in tiff format [12] using the compression lossless method. Third, histogram stretching (HS) was employed to enhance image contrast. Suppose the kth input and output of HS is xk and yk. HS can be formulated as
(1)yk=xk−xminkxmaxk−xmink,where xmink and xmaxk stand for the minimum and maximum grayscale values in the input xk.Fourth, cropping was done in order to eliminate the check-up bed at the bottom, the subject’s two arms at bilateral sides, the rulers at the bottom and right side, and information (hospital, scanning protocol, subject’s information, image head information, and labeling) at four corners.Lastly, downscaling was performed to reduce each image to the size of256×256. Figure 1 displays the diagram of our preprocessing procedure. Figures 2(a) and 2(b) shows two preprocessed examples of TOF and HC, respectively.Figure 1
Diagram of preprocessing.Figure 2
Illustration of our dataset.
(a)
TOF(b)
HC


### 3.2. Stochastic Pooling

Pooling is an essential operation in standard convolutional neural networks (CNNs) [13]. Two types of pooling exist. One is max pooling (MP), and the other is average pooling (AP). The objective of pooling is to down-sample an input image or feature map (FM), reducing their dimensionality (width or height) and allowing for some assumption about the features to be made in each block.Suppose we have an input image or FM, which can be split intoO1×O2 blocks, where every block has the extent of Q1×Q2. Currently, let us fix on the block Bo1,o2 at o1th row and o2th column as shown as the red rectangle in Figure 3.
(2)Bo1,o2=bx,y,x=1,⋯,Q1,y=1,⋯,Q2,where 1≤o1≤O1,1≤o2≤O2, bx,y means the pixel value at coordinate x,y.Figure 3
A diagram of block-wise pooling.The strided convolution (SC) goes over the input activation map with the strides that equals the size of the blockQ1,Q2. The output of SC is
(3)Bo1,o2SC=b1,1.Thel2-norm pooling (L2P), average pooling (AP) [14], and max pooling (MP) [15] produce the l2-norm, average, and maximum values within the block Bm1,m2, respectively. Their formula can be written as below:
(4)Bo1,o2L2P=∑x=1Q1∑y=1Q2b2x,yQ1×Q2,(5)Bo1,o2AP=1Q1×Q2∑x=1Q1∑y=1Q2bx,y,(6)Bo1,o2MP=maxx=1Q1maxy=1Q2bx,y.Nevertheless, the AP outputs the average, downscaling the greatest value, where the important features may lie. In contrast, MP stores the greatest value but deteriorates the overfitting obstacle. In order to solve the above concerns, stochastic pooling (SP) [15] is introduced to provide a resolution to the drawbacks of AP and MP. SP is a four-step process.Step 1.
It produces the probability map (PM)VPM for each pixel in the block Bo1,o2.
(7)VPMx,y=bx,y∑x=1Q1∑y=1Q2bx,y,s.t.∑x=1Q1∑y=1Q2VPMx,y=1,where VPMx,y stands for the PM value at pixel x,y.Step 2.
It creates a random location vector (RLV)r⟶=xr,yr that takes the discrete probability distribution (DPD) as
(8)Pr⟶=1,1=VPM1,1Pr⟶=2,1=VPM2,1⋯Pr⟶=Q1,1=VPMQ1,1Pr⟶=1,2=VPM1,2⋯Pr⟶=1,Q2=VPM1,Q2,Prr⟶=2,2=VPM2,2⋯Pr⟶=2,Q2=VPM2,Q2,⋯⋯⋯Pr⟶=Q1,2=VPMQ1,2⋯Pr⟶=Q1,Q2=VPMQ1,Q2,where P represents the probability. Shortly speaking, Pr⟶=x,y=VPMx,y, ∀1≤x≤Q1&1≤y≤Q2 or r⟶~VPM, namely, the distribution of RLV r⟶ has the DPD as VPM.Step 3.
A sample location vectorr⟶0 is drawn from the RLV r⟶, and we have
(9)r⟶0=xr0,yr0.Step 4.
SP outputs the pixel at the locationr0⟶, namely,
(10)Bo1,o2SP=bxr0,yr0.Algorithm 1: Pseudocode of SP.
Input: blockBo1,o2.Step 1: produce the PM for each pixelVPM. See Equation (7).Step 2: create a RLVr⟶=xr,yr. See Equation (8).Step 3: draw a sample location vectorr⟶0 from the RLV r⟶. See Equation (9).Step 4: output pixel at locationr⟶0. See Equation (10)Output matrixBo1,o2SP.Figure4 shows a realistic example of four different pooling methods. Algorithm 1 presents the pseudocode of SP. Take the 3×3 block B1,1 (The red rectangle in Figure 4) as an example, L2P generates the output as 6.98. AP and MP present 5.99 and 9.9, respectively. Meanwhile, SP first generates PM matrix
(11)VPM=0.180.180.030.050.10.010.090.180.18,and a sample location vector is drawn as r⟶0=1,2. Therefore, the output of SP is B1,1SP=b1,2=9.8.Figure 4
Comparison of four different pooling methods.


### 3.3. Structural Optimization

How to obtain the best network structure [16]? We try to design nine different configurations in this study. Their hyperparameters of structures are listed in Table 2. Two hyperparameters are considered in this study: (i) the number of Conv layers and (ii) the number of fully connected layers (FCLs). Those two types of layers are common layers in standard CNN, so we will not introduce them due to the page limit.Table 2
Structures of nine customized neural networks.
ConfigurationNo. of Conv layersNo. of FCLsI21II22III23IV31V32VI33VII41VIII42IX43Bold means the best.In the following experiment, we will observe that configuration V, a five-layer customized neural network, gives the best performance. Here, we briefly give its detailed structure in Table3. The input is of size 256×256×1. The first Conv layer (Conv_1) is associated with the batch normalization (BN) layer and rectified linear unit (ReLU) activation. The parameters of Conv_1 are 32 kernels with sizes of 3×3 and stride of 2. Afterward, the first SP (SP_1) reduce the FM from 128×128×32 to 64×64×32.Table 3
Detailed structure of network of configuration V.
LayerParametersFMInput256×256×1Conv_1 (BN-ReLU)32,3×3/2128×128×32SP_164×64×32Conv_2 (BN-ReLU)64,3×3/232×32×64SP_216×16×64Conv_3 (BN-ReLU)128,3×316×16×128SP_38×8×128Flatten8192FCL_1100×8192, 100×1100FCL_22×100,2×12OutputAfter three Conv layers and three SP layers, the size of FM is8×8×128. It is then flattened to a vector of 8192 neurons. With two FCLs of 100 and 2 hidden neurons, the neural network finally outputs whether TOF or HC. All in all, our model is termed structurally optimized stochastic pooling convolutional neural network (SOSPCNN). The FM plot is portrayed in Figure 5.Figure 5
FM plot.


### 3.4. Multiple-Way Data Augmentation

The relatively small dataset (40+40 = 80 images) may bring the overfitting problem. To avoid overfitting, data augmentation (DA) [17] is a powerful tool because it can generate synthetic images on the training set [18]. Zhu (2021) [19] presented an 18-way DA method and proved this 18-way DA works better than the traditional DA approach. Its diagram is shown in Figure 6. The difference of DA and MDA is that (i) MDA uses a combination of different DA methods on training set; (ii) MDA is modular design. The users are easy to add or remove particular DA methods from a MDA.Figure 6
Diagram of multiple-way data augmentation.Suppose we have the raw training imagerk, where k represents the image index. First, P1 different DA methods displayed in Figure 6 are applied to rk. Let Zp,p=1,⋯,P1 be each DA operation, we get P1 augmented datasets on raw image rk as
(12)Zprk,p=1,⋯,P1.LetP2 stands for the size of generated new images for each DA method, thus,
(13)Zprk=P2.Second, horizontal mirrored image is produced by(14)rhk=β1rk,where β1 means horizontal mirror function.Third, allP1 different DA methods are carried out on the mirrored image rhk and produce P1 new datasets as
(15)Zprhk,p=1,⋯,P1,Zprhk=P2,p=1,⋯,P1.Fourth, the raw imagerk, the mirrored image rhk, all P1-way results of raw image Zprk, and all P1-way DA results of horizontal mirrored image Zprhk are combined. The final generated dataset from rk is defined as Gk:
(16)rk↦Gk=β2rkrhkZ1rk⏟P2Z1rhk⏟P2⋯⋯ZP1rk⏟P2ZP1rhk⏟P2,where β2 stands for the concatenation function. Let augmentation factor be P3, which means the number of images in Gk, we obtain
(17)P3=Gkrk=1+P1×P2×21=2×P1×P2+2.Algorithm2 recaps the pseudocode of this 18-way DA. We set P1=9 to achieve an 18-way DA. We also set P2=30, thus P3=542, indicating each raw training image will generate 542 images, which include the raw image rk itself.Algorithm 2: Pseudocode of proposed 18-way DA on kth training image
Input: import raw preprocessedkth training image rk.P1 geometric or photometric or noise-injection DA transforms Zp are utilized on rk.Step 1: we obtainZprk,p=1,⋯,P1. See Equation (12)Each enhanced dataset containsP2 new images. See Equation (13).Step 2: a horizontal mirror image is produced asrhk=β1rk. See Equation (14).Step 3:M1-way data augmentation methods are carried out on rhk,We obtainZprhk,p=1,⋯,P1. See Equation (15).Step 4:rk, rhk, Zprk,p=1,⋯,P1, and Zprhk,p=1,⋯,P1 are merged via β2. See Equation (16).Output A new datasetGk is produced. The number of images is P3=2×P1×P2+2. See Equation (17).


### 3.5. Implementation and Grad-CAM

Q-fold crossvalidation [20] is employed. The whole dataset is divided into Q folds (see Figure 7). At qth trial, 1≤q≤Q, the qth fold is picked up as the test, and the rest Q−1 folds: 1,⋯,q−1,q+1,⋯,Q are chosen as training set [21]. In this study, we set Q=10, namely, a 10-fold cross validation. Furthermore, we run the 10-fold crossvalidation 10 times, i.e., 10×10-fold crossvalidation.Figure 7
Q-fold crossvalidation.Gradient-weighted class activation mapping (Grad-CAM) [22] is employed to explain how our model makes its decision in classification. Grad-CAM utilizes the gradient of the classification score with respect to the convolutional features determined by the network to understand which parts of the image are most important for classification. Grad-CAM is a generalization of the class activation mapping (CAM) method [23] to a broader range of CNN models since the original CAM relies on a fully convolutional neural network structure. The output of SP_3 (see Table 3) is used as the feature layer for Grad-CAM.Mathematically, suppose our classification network is with outputyc, standing for the score for class c. We would like to compute the Grad-CAM map for a layer with k feature maps Ai,jk, where i,j stands for the indexes of pixels. We can obtain the neural importance weight as
(18)αkc=1N∑i∑j∂yc∂Ai,jk,where N stands for the total number of pixels in the feature map. The Grad-CAM is a weighted combination of the feature maps with a ReLU as
(19)M=ReLU∑kαkcAk.The Grad-CAM mapM is then upsampled to the size of input data.


### 3.6. Measures

The confusion matrix of 10 runs of 10-fold crossvalidation is supposed to be(20)G=g1,1g1,2g2,1g2,2=TPFNFPTN.NoteFN=FP=0 for a perfect classification. The meaning of P, N, TP, FP, TN, and FN are itemized in Table 4.Table 4
Meanings in measures.
AbbreviationFull formSymbolMeaningPPositiveTOFNNegativeHCTPTrue positiveg1,1TOF images are classified correctly.FPFalse positiveg2,1HC images are wrongly classified as TOF.TNTrue negativeg2,2HC images are classified correctly.FNFalse negativeg1,2TOF images are wrongly classified as HC.Nine measures are used: sensitivity, specificity, precision, accuracy, F1 score, Matthews correlation coefficient (MCC), Fowlkes–Mallows index (FMI), receiver operating characteristic (ROC), and area under the curve (AUC).The first four measures are defined as(21)Sen=g1,1g1,1+g1,2Spc=g2,2g2,2+g2,1,Prc=g1,1g1,1+g2,1Acc=g1,1+g2,2g1,1+g2,2+g1,2+g2,1.F1, MCC [24], and FMI [25] are defined as
(22)F1=2×Sen×PrcSen+Prc=2×g1,12×g1,1+g1,2+g2,1,MCC=g1,1×g2,2−g2,1×g1,2g1,1+g2,1×g1,1+g1,2×g2,2+g2,1×g2,2+g1,2,FMI=Sen×Prc=g1,1g1,1+g1,2×g1,1g1,1+g2,1.The above measures are calculated in the mean and standard deviation (MSD) format. Furthermore, ROC is a curve to measure a binary classifier with varying discrimination thresholds [26]. The ROC curve is created by plotting the sensitivity against 1-specificity. The AUC is calculated based on the ROC curve [27].


## 3.1. Preprocessing

Table1 lists the abbreviation list for the ease of reading. A five-step preprocessing was carried out on all the images to select the important slices, save storage, enhance contrast, remove unnecessary image regions, and reduce the image resolution.Table 1
Abbreviation list.
AbbreviationMeaningAPAverage poolingAUCArea under the curveBNNormalizationCNNConvolutional neural networkCTComputed tomographyDAData augmentationDPDDiscrete probability distributionFCLFully connected layerFMIFowlkes–Mallows indexGrad-CAMGradient-weighted class activation mappingGUIGraphical user interfaceHSHistogram stretchingL2Pl2-norm poolingMCCMatthews correlation coefficientMPMax-poolingMRIMagnetic resonance imagingMSDMean and standard deviationPMProbability mapReLURectified linear unitRLVRandom location vectorROCReceiver operating characteristicSCStrided convolutionSPStochastic poolingTOFTetralogy of FallotFirst, four slices were chosen by radiologists using a slice-level selection method. For TOF patients, the slices showing the largest size and number of lesions were selected. For healthy control subjects, any level of the image can be selected. Now, we have in total 40 TOF images and 40 HC images.Second, all the images are converted to grayscale images and stored in tiff format [12] using the compression lossless method. Third, histogram stretching (HS) was employed to enhance image contrast. Suppose the kth input and output of HS is xk and yk. HS can be formulated as
(1)yk=xk−xminkxmaxk−xmink,where xmink and xmaxk stand for the minimum and maximum grayscale values in the input xk.Fourth, cropping was done in order to eliminate the check-up bed at the bottom, the subject’s two arms at bilateral sides, the rulers at the bottom and right side, and information (hospital, scanning protocol, subject’s information, image head information, and labeling) at four corners.Lastly, downscaling was performed to reduce each image to the size of256×256. Figure 1 displays the diagram of our preprocessing procedure. Figures 2(a) and 2(b) shows two preprocessed examples of TOF and HC, respectively.Figure 1
Diagram of preprocessing.Figure 2
Illustration of our dataset.
(a)
TOF(b)
HC


## 3.2. Stochastic Pooling

Pooling is an essential operation in standard convolutional neural networks (CNNs) [13]. Two types of pooling exist. One is max pooling (MP), and the other is average pooling (AP). The objective of pooling is to down-sample an input image or feature map (FM), reducing their dimensionality (width or height) and allowing for some assumption about the features to be made in each block.Suppose we have an input image or FM, which can be split intoO1×O2 blocks, where every block has the extent of Q1×Q2. Currently, let us fix on the block Bo1,o2 at o1th row and o2th column as shown as the red rectangle in Figure 3.
(2)Bo1,o2=bx,y,x=1,⋯,Q1,y=1,⋯,Q2,where 1≤o1≤O1,1≤o2≤O2, bx,y means the pixel value at coordinate x,y.Figure 3
A diagram of block-wise pooling.The strided convolution (SC) goes over the input activation map with the strides that equals the size of the blockQ1,Q2. The output of SC is
(3)Bo1,o2SC=b1,1.Thel2-norm pooling (L2P), average pooling (AP) [14], and max pooling (MP) [15] produce the l2-norm, average, and maximum values within the block Bm1,m2, respectively. Their formula can be written as below:
(4)Bo1,o2L2P=∑x=1Q1∑y=1Q2b2x,yQ1×Q2,(5)Bo1,o2AP=1Q1×Q2∑x=1Q1∑y=1Q2bx,y,(6)Bo1,o2MP=maxx=1Q1maxy=1Q2bx,y.Nevertheless, the AP outputs the average, downscaling the greatest value, where the important features may lie. In contrast, MP stores the greatest value but deteriorates the overfitting obstacle. In order to solve the above concerns, stochastic pooling (SP) [15] is introduced to provide a resolution to the drawbacks of AP and MP. SP is a four-step process.Step 1.
It produces the probability map (PM)VPM for each pixel in the block Bo1,o2.
(7)VPMx,y=bx,y∑x=1Q1∑y=1Q2bx,y,s.t.∑x=1Q1∑y=1Q2VPMx,y=1,where VPMx,y stands for the PM value at pixel x,y.Step 2.
It creates a random location vector (RLV)r⟶=xr,yr that takes the discrete probability distribution (DPD) as
(8)Pr⟶=1,1=VPM1,1Pr⟶=2,1=VPM2,1⋯Pr⟶=Q1,1=VPMQ1,1Pr⟶=1,2=VPM1,2⋯Pr⟶=1,Q2=VPM1,Q2,Prr⟶=2,2=VPM2,2⋯Pr⟶=2,Q2=VPM2,Q2,⋯⋯⋯Pr⟶=Q1,2=VPMQ1,2⋯Pr⟶=Q1,Q2=VPMQ1,Q2,where P represents the probability. Shortly speaking, Pr⟶=x,y=VPMx,y, ∀1≤x≤Q1&1≤y≤Q2 or r⟶~VPM, namely, the distribution of RLV r⟶ has the DPD as VPM.Step 3.
A sample location vectorr⟶0 is drawn from the RLV r⟶, and we have
(9)r⟶0=xr0,yr0.Step 4.
SP outputs the pixel at the locationr0⟶, namely,
(10)Bo1,o2SP=bxr0,yr0.Algorithm 1: Pseudocode of SP.
Input: blockBo1,o2.Step 1: produce the PM for each pixelVPM. See Equation (7).Step 2: create a RLVr⟶=xr,yr. See Equation (8).Step 3: draw a sample location vectorr⟶0 from the RLV r⟶. See Equation (9).Step 4: output pixel at locationr⟶0. See Equation (10)Output matrixBo1,o2SP.Figure4 shows a realistic example of four different pooling methods. Algorithm 1 presents the pseudocode of SP. Take the 3×3 block B1,1 (The red rectangle in Figure 4) as an example, L2P generates the output as 6.98. AP and MP present 5.99 and 9.9, respectively. Meanwhile, SP first generates PM matrix
(11)VPM=0.180.180.030.050.10.010.090.180.18,and a sample location vector is drawn as r⟶0=1,2. Therefore, the output of SP is B1,1SP=b1,2=9.8.Figure 4
Comparison of four different pooling methods.


## 3.3. Structural Optimization

How to obtain the best network structure [16]? We try to design nine different configurations in this study. Their hyperparameters of structures are listed in Table 2. Two hyperparameters are considered in this study: (i) the number of Conv layers and (ii) the number of fully connected layers (FCLs). Those two types of layers are common layers in standard CNN, so we will not introduce them due to the page limit.Table 2
Structures of nine customized neural networks.
ConfigurationNo. of Conv layersNo. of FCLsI21II22III23IV31V32VI33VII41VIII42IX43Bold means the best.In the following experiment, we will observe that configuration V, a five-layer customized neural network, gives the best performance. Here, we briefly give its detailed structure in Table3. The input is of size 256×256×1. The first Conv layer (Conv_1) is associated with the batch normalization (BN) layer and rectified linear unit (ReLU) activation. The parameters of Conv_1 are 32 kernels with sizes of 3×3 and stride of 2. Afterward, the first SP (SP_1) reduce the FM from 128×128×32 to 64×64×32.Table 3
Detailed structure of network of configuration V.
LayerParametersFMInput256×256×1Conv_1 (BN-ReLU)32,3×3/2128×128×32SP_164×64×32Conv_2 (BN-ReLU)64,3×3/232×32×64SP_216×16×64Conv_3 (BN-ReLU)128,3×316×16×128SP_38×8×128Flatten8192FCL_1100×8192, 100×1100FCL_22×100,2×12OutputAfter three Conv layers and three SP layers, the size of FM is8×8×128. It is then flattened to a vector of 8192 neurons. With two FCLs of 100 and 2 hidden neurons, the neural network finally outputs whether TOF or HC. All in all, our model is termed structurally optimized stochastic pooling convolutional neural network (SOSPCNN). The FM plot is portrayed in Figure 5.Figure 5
FM plot.


## 3.4. Multiple-Way Data Augmentation

The relatively small dataset (40+40 = 80 images) may bring the overfitting problem. To avoid overfitting, data augmentation (DA) [17] is a powerful tool because it can generate synthetic images on the training set [18]. Zhu (2021) [19] presented an 18-way DA method and proved this 18-way DA works better than the traditional DA approach. Its diagram is shown in Figure 6. The difference of DA and MDA is that (i) MDA uses a combination of different DA methods on training set; (ii) MDA is modular design. The users are easy to add or remove particular DA methods from a MDA.Figure 6
Diagram of multiple-way data augmentation.Suppose we have the raw training imagerk, where k represents the image index. First, P1 different DA methods displayed in Figure 6 are applied to rk. Let Zp,p=1,⋯,P1 be each DA operation, we get P1 augmented datasets on raw image rk as
(12)Zprk,p=1,⋯,P1.LetP2 stands for the size of generated new images for each DA method, thus,
(13)Zprk=P2.Second, horizontal mirrored image is produced by(14)rhk=β1rk,where β1 means horizontal mirror function.Third, allP1 different DA methods are carried out on the mirrored image rhk and produce P1 new datasets as
(15)Zprhk,p=1,⋯,P1,Zprhk=P2,p=1,⋯,P1.Fourth, the raw imagerk, the mirrored image rhk, all P1-way results of raw image Zprk, and all P1-way DA results of horizontal mirrored image Zprhk are combined. The final generated dataset from rk is defined as Gk:
(16)rk↦Gk=β2rkrhkZ1rk⏟P2Z1rhk⏟P2⋯⋯ZP1rk⏟P2ZP1rhk⏟P2,where β2 stands for the concatenation function. Let augmentation factor be P3, which means the number of images in Gk, we obtain
(17)P3=Gkrk=1+P1×P2×21=2×P1×P2+2.Algorithm2 recaps the pseudocode of this 18-way DA. We set P1=9 to achieve an 18-way DA. We also set P2=30, thus P3=542, indicating each raw training image will generate 542 images, which include the raw image rk itself.Algorithm 2: Pseudocode of proposed 18-way DA on kth training image
Input: import raw preprocessedkth training image rk.P1 geometric or photometric or noise-injection DA transforms Zp are utilized on rk.Step 1: we obtainZprk,p=1,⋯,P1. See Equation (12)Each enhanced dataset containsP2 new images. See Equation (13).Step 2: a horizontal mirror image is produced asrhk=β1rk. See Equation (14).Step 3:M1-way data augmentation methods are carried out on rhk,We obtainZprhk,p=1,⋯,P1. See Equation (15).Step 4:rk, rhk, Zprk,p=1,⋯,P1, and Zprhk,p=1,⋯,P1 are merged via β2. See Equation (16).Output A new datasetGk is produced. The number of images is P3=2×P1×P2+2. See Equation (17).


## 3.5. Implementation and Grad-CAM

Q-fold crossvalidation [20] is employed. The whole dataset is divided into Q folds (see Figure 7). At qth trial, 1≤q≤Q, the qth fold is picked up as the test, and the rest Q−1 folds: 1,⋯,q−1,q+1,⋯,Q are chosen as training set [21]. In this study, we set Q=10, namely, a 10-fold cross validation. Furthermore, we run the 10-fold crossvalidation 10 times, i.e., 10×10-fold crossvalidation.Figure 7
Q-fold crossvalidation.Gradient-weighted class activation mapping (Grad-CAM) [22] is employed to explain how our model makes its decision in classification. Grad-CAM utilizes the gradient of the classification score with respect to the convolutional features determined by the network to understand which parts of the image are most important for classification. Grad-CAM is a generalization of the class activation mapping (CAM) method [23] to a broader range of CNN models since the original CAM relies on a fully convolutional neural network structure. The output of SP_3 (see Table 3) is used as the feature layer for Grad-CAM.Mathematically, suppose our classification network is with outputyc, standing for the score for class c. We would like to compute the Grad-CAM map for a layer with k feature maps Ai,jk, where i,j stands for the indexes of pixels. We can obtain the neural importance weight as
(18)αkc=1N∑i∑j∂yc∂Ai,jk,where N stands for the total number of pixels in the feature map. The Grad-CAM is a weighted combination of the feature maps with a ReLU as
(19)M=ReLU∑kαkcAk.The Grad-CAM mapM is then upsampled to the size of input data.


## 3.6. Measures

The confusion matrix of 10 runs of 10-fold crossvalidation is supposed to be(20)G=g1,1g1,2g2,1g2,2=TPFNFPTN.NoteFN=FP=0 for a perfect classification. The meaning of P, N, TP, FP, TN, and FN are itemized in Table 4.Table 4
Meanings in measures.
AbbreviationFull formSymbolMeaningPPositiveTOFNNegativeHCTPTrue positiveg1,1TOF images are classified correctly.FPFalse positiveg2,1HC images are wrongly classified as TOF.TNTrue negativeg2,2HC images are classified correctly.FNFalse negativeg1,2TOF images are wrongly classified as HC.Nine measures are used: sensitivity, specificity, precision, accuracy, F1 score, Matthews correlation coefficient (MCC), Fowlkes–Mallows index (FMI), receiver operating characteristic (ROC), and area under the curve (AUC).The first four measures are defined as(21)Sen=g1,1g1,1+g1,2Spc=g2,2g2,2+g2,1,Prc=g1,1g1,1+g2,1Acc=g1,1+g2,2g1,1+g2,2+g1,2+g2,1.F1, MCC [24], and FMI [25] are defined as
(22)F1=2×Sen×PrcSen+Prc=2×g1,12×g1,1+g1,2+g2,1,MCC=g1,1×g2,2−g2,1×g1,2g1,1+g2,1×g1,1+g1,2×g2,2+g2,1×g2,2+g1,2,FMI=Sen×Prc=g1,1g1,1+g1,2×g1,1g1,1+g2,1.The above measures are calculated in the mean and standard deviation (MSD) format. Furthermore, ROC is a curve to measure a binary classifier with varying discrimination thresholds [26]. The ROC curve is created by plotting the sensitivity against 1-specificity. The AUC is calculated based on the ROC curve [27].


## 4. Experimental Results

### 4.1. Statistical Analysis

The result of the SOSPCNN model using configuration V is itemized in Table5. The model arrives at a performance with a sensitivity of 92.25±2.19, a specificity of 92.75±2.49, a precision of 92.79±2.29, an accuracy of 92.50±1.18, an F1 score of 92.48±1.17, an MCC of 85.06±2.38, and an FMI of 92.50±1.17.Table 5
Statistical analysis of SOSPCNN model.
RunSenSpcPrcAccF1MCCFMI195.0092.5092.6893.7593.8387.5393.83292.5090.0090.2491.2591.3682.5391.36395.0092.5092.6893.7593.8387.5393.83490.0092.5092.3191.2591.1482.5391.15590.0095.0094.7492.5092.3185.1192.34690.0097.5097.3093.7593.5187.7593.58792.5095.0094.8793.7593.6787.5393.68892.5090.0090.2491.2591.3682.5391.36990.0092.5092.3191.2591.1482.5391.151095.0090.0090.4892.5092.6885.1192.71MSD92.25±2.1992.75±2.4992.79±2.2992.50±1.1892.48±1.1785.06±2.3892.50±1.17Figure8 shows the confusion matrix of 10×10-fold crossvalidation, where we can see the TP=369, FN=31, TN=371, and FP=29, indicating 31 TOF are wrongly classified as HC while 29 HC are misclassified to TOF. Hence, the sensitivity is 369/369+31=92.25%, and specificity is 371/29+371=92.75%.Figure 8
Confusion matrix of10×10-fold crossvalidation (Here, classes 1 and 2 stand for ToF and HC, respectively).


### 4.2. Configuration Comparison

We compare nine configurations (see Table2). The validation is the same as previous experiment. Due to the page limit, the detailed statistical analysis is not shown. The ROC and AUC values are displayed in Figure 9. The AUC values of nine networks with different configurations are: 0.9502, 0.9511, 0.9504, 0.9532, 0.9587, 0.9577, 0.9360, 0.9419, and 0.9389 (as shown in Figure 10). We can observe from Figure 10 that the best network is with configuration V, whose structure is shown in Table 3.Figure 9
Comparison of nine configurations.
(a)
Configuration I(b)
Configuration II(c)
Configuration III(d)
Configuration IV(e)
Configuration V(f)
Configuration VI(g)
Configuration VII(h)
Configuration VIII(i)
Configuration IXFigure 10
Bar plot of AUC against nine configurations.


### 4.3. Effect of Multiple-Way Data Augmentation

Figure11 shows the multiple-way DA results if we take Figure 2(a) as the raw training examples. Due to the page limit, the multiple-way DA results on the horizontally mirrored image are not displayed. As we can see from Figure 11, multiple-way DA increases the diversity of the training images.Figure 11
Illustration of multiple-way data augmentation.
(a)
Gaussian noise(b)
Salt-and-pepper noise(c)
Speckle noise(d)
Horizontal shear(e)
Vertical shear(f)
Rotation(g)
Gamma correction(h)
Random translation(i)
ScalingIf we remove the multiple-way data augmentation from our model, the performances are decreased, as shown in Table6, where MSD stands for mean and standard deviation. Comparing Table 5 with Table 6, we can observe multiple-way DA is efficient in improving the classification performance. The reason is that it helps our model resist overfitting by enhancing the diversity of the training set.Table 6
Statistical analysis without multiple-way data augmentation.
RunSenSpcPrcAccF1MCCFMI187.5090.0089.7488.7588.6177.5288.61290.0087.5087.8088.7588.8977.5288.90387.5087.5087.5087.5087.5075.0087.50490.0087.5087.8088.7588.8977.5288.90585.0090.0089.4787.5087.1875.0987.21685.0090.0089.4787.5087.1875.0987.21782.5092.5091.6787.5086.8475.3886.96887.5087.5087.5087.5087.5075.0087.50985.0090.0089.4787.5087.1875.0987.211082.5092.5091.6787.5086.8475.3886.96MSD86.25 ± 2.7089.50 ± 1.9789.21 ± 1.5887.87 ± 0.6087.66 ± 0.8275.86 ± 1.1687.70 ± 0.79


### 4.4. Explainability

Figure12 shows the manual delineation and the heat map of Figure 2(a) via Grad-CAM described in Section 3.5. The manual delineation showed the radiologist make decisions on “TOF” diagnosis based on all the areas of the abnormal heart, while the heat map shows the proposed SOSPCNN model also puts more focus on the heart region other than the surrounding tissues and background areas.Figure 12
Heatmap of one TOF image.
(a)
Manually delineated(b)
Heatmap


### 4.5. Comparison with State-of-the-Art Approaches

We compare the proposed SOSPCNN model with three other approaches: MC [8], 3DCNN [9], and VCCNN [10]. The results are shown in Table 7. Note that some comparison methods are not suitable for our dataset, so we modify them to adapt to our dataset.Table 7
Comparison with state-of-the-art approaches.
ApproachSenSpcPrcAccF1MCCFMIMC [8]86.25±3.5880.75±3.5581.88±2.3283.50±0.7983.92±0.9767.25±1.5684.00±0.983DCNN [9]91.00±3.1689.50±3.2989.77±2.7090.25±1.4290.32±1.4280.63±2.7990.35±1.41VCCNN [10]90.75±1.6990.00±2.3690.14±1.9590.38±0.8490.41±0.7880.80±1.6490.43±0.76SOSPCNN (ours)92.25±2.1992.75±2.4992.79±2.2992.50±1.1892.48±1.1785.06±2.3892.50±1.17The error bar comparison is drawn in Figure13, which clearly shows that the proposed SOSPCNN outperforms all three comparative approaches. The reason is three folds: (i) we use stochastic pooling to replace traditional max-pooling; (ii) we use structural optimization to determine the optimal structure of our SOSPCNN model; (iii) multiple-way DA is included to increase the diversity of training images. In the future, more advanced techniques [28–30] will be tested and integrated into our model.Figure 13
Error bar comparison.


### 4.6. Desktop and Web Apps

MATLAB app designer is used to create a professional application for both desktop and web. The input to this web app is any cardiovascular CT image, and the aforementioned SOSPCNN model is included in our app. Figure14(a) displays the graphical user interface (GUI) of the standalone desktop app. The users can upload their custom images, and the software can show the results by turning the knob into the correct texts: TOF, HC, or none.Figure 14
GUI of developed apps.
(a)
Desktop app(b)
Web appFigure14(b) shows the GUI of the web app that is accessed through a “Google Chrome” web browser. The web app is based on a client-server modeled structure [31], i.e., the user is provided services through an off-site server hosted by a third-party cloud service, Microsoft Azure in our study. Our developed online web app can assist hospital clinicians in making decisions remotely and effectively.


## 4.1. Statistical Analysis

The result of the SOSPCNN model using configuration V is itemized in Table5. The model arrives at a performance with a sensitivity of 92.25±2.19, a specificity of 92.75±2.49, a precision of 92.79±2.29, an accuracy of 92.50±1.18, an F1 score of 92.48±1.17, an MCC of 85.06±2.38, and an FMI of 92.50±1.17.Table 5
Statistical analysis of SOSPCNN model.
RunSenSpcPrcAccF1MCCFMI195.0092.5092.6893.7593.8387.5393.83292.5090.0090.2491.2591.3682.5391.36395.0092.5092.6893.7593.8387.5393.83490.0092.5092.3191.2591.1482.5391.15590.0095.0094.7492.5092.3185.1192.34690.0097.5097.3093.7593.5187.7593.58792.5095.0094.8793.7593.6787.5393.68892.5090.0090.2491.2591.3682.5391.36990.0092.5092.3191.2591.1482.5391.151095.0090.0090.4892.5092.6885.1192.71MSD92.25±2.1992.75±2.4992.79±2.2992.50±1.1892.48±1.1785.06±2.3892.50±1.17Figure8 shows the confusion matrix of 10×10-fold crossvalidation, where we can see the TP=369, FN=31, TN=371, and FP=29, indicating 31 TOF are wrongly classified as HC while 29 HC are misclassified to TOF. Hence, the sensitivity is 369/369+31=92.25%, and specificity is 371/29+371=92.75%.Figure 8
Confusion matrix of10×10-fold crossvalidation (Here, classes 1 and 2 stand for ToF and HC, respectively).


## 4.2. Configuration Comparison

We compare nine configurations (see Table2). The validation is the same as previous experiment. Due to the page limit, the detailed statistical analysis is not shown. The ROC and AUC values are displayed in Figure 9. The AUC values of nine networks with different configurations are: 0.9502, 0.9511, 0.9504, 0.9532, 0.9587, 0.9577, 0.9360, 0.9419, and 0.9389 (as shown in Figure 10). We can observe from Figure 10 that the best network is with configuration V, whose structure is shown in Table 3.Figure 9
Comparison of nine configurations.
(a)
Configuration I(b)
Configuration II(c)
Configuration III(d)
Configuration IV(e)
Configuration V(f)
Configuration VI(g)
Configuration VII(h)
Configuration VIII(i)
Configuration IXFigure 10
Bar plot of AUC against nine configurations.


## 4.3. Effect of Multiple-Way Data Augmentation

Figure11 shows the multiple-way DA results if we take Figure 2(a) as the raw training examples. Due to the page limit, the multiple-way DA results on the horizontally mirrored image are not displayed. As we can see from Figure 11, multiple-way DA increases the diversity of the training images.Figure 11
Illustration of multiple-way data augmentation.
(a)
Gaussian noise(b)
Salt-and-pepper noise(c)
Speckle noise(d)
Horizontal shear(e)
Vertical shear(f)
Rotation(g)
Gamma correction(h)
Random translation(i)
ScalingIf we remove the multiple-way data augmentation from our model, the performances are decreased, as shown in Table6, where MSD stands for mean and standard deviation. Comparing Table 5 with Table 6, we can observe multiple-way DA is efficient in improving the classification performance. The reason is that it helps our model resist overfitting by enhancing the diversity of the training set.Table 6
Statistical analysis without multiple-way data augmentation.
RunSenSpcPrcAccF1MCCFMI187.5090.0089.7488.7588.6177.5288.61290.0087.5087.8088.7588.8977.5288.90387.5087.5087.5087.5087.5075.0087.50490.0087.5087.8088.7588.8977.5288.90585.0090.0089.4787.5087.1875.0987.21685.0090.0089.4787.5087.1875.0987.21782.5092.5091.6787.5086.8475.3886.96887.5087.5087.5087.5087.5075.0087.50985.0090.0089.4787.5087.1875.0987.211082.5092.5091.6787.5086.8475.3886.96MSD86.25 ± 2.7089.50 ± 1.9789.21 ± 1.5887.87 ± 0.6087.66 ± 0.8275.86 ± 1.1687.70 ± 0.79


## 4.4. Explainability

Figure12 shows the manual delineation and the heat map of Figure 2(a) via Grad-CAM described in Section 3.5. The manual delineation showed the radiologist make decisions on “TOF” diagnosis based on all the areas of the abnormal heart, while the heat map shows the proposed SOSPCNN model also puts more focus on the heart region other than the surrounding tissues and background areas.Figure 12
Heatmap of one TOF image.
(a)
Manually delineated(b)
Heatmap


## 4.5. Comparison with State-of-the-Art Approaches

We compare the proposed SOSPCNN model with three other approaches: MC [8], 3DCNN [9], and VCCNN [10]. The results are shown in Table 7. Note that some comparison methods are not suitable for our dataset, so we modify them to adapt to our dataset.Table 7
Comparison with state-of-the-art approaches.
ApproachSenSpcPrcAccF1MCCFMIMC [8]86.25±3.5880.75±3.5581.88±2.3283.50±0.7983.92±0.9767.25±1.5684.00±0.983DCNN [9]91.00±3.1689.50±3.2989.77±2.7090.25±1.4290.32±1.4280.63±2.7990.35±1.41VCCNN [10]90.75±1.6990.00±2.3690.14±1.9590.38±0.8490.41±0.7880.80±1.6490.43±0.76SOSPCNN (ours)92.25±2.1992.75±2.4992.79±2.2992.50±1.1892.48±1.1785.06±2.3892.50±1.17The error bar comparison is drawn in Figure13, which clearly shows that the proposed SOSPCNN outperforms all three comparative approaches. The reason is three folds: (i) we use stochastic pooling to replace traditional max-pooling; (ii) we use structural optimization to determine the optimal structure of our SOSPCNN model; (iii) multiple-way DA is included to increase the diversity of training images. In the future, more advanced techniques [28–30] will be tested and integrated into our model.Figure 13
Error bar comparison.


## 4.6. Desktop and Web Apps

MATLAB app designer is used to create a professional application for both desktop and web. The input to this web app is any cardiovascular CT image, and the aforementioned SOSPCNN model is included in our app. Figure14(a) displays the graphical user interface (GUI) of the standalone desktop app. The users can upload their custom images, and the software can show the results by turning the knob into the correct texts: TOF, HC, or none.Figure 14
GUI of developed apps.
(a)
Desktop app(b)
Web appFigure14(b) shows the GUI of the web app that is accessed through a “Google Chrome” web browser. The web app is based on a client-server modeled structure [31], i.e., the user is provided services through an off-site server hosted by a third-party cloud service, Microsoft Azure in our study. Our developed online web app can assist hospital clinicians in making decisions remotely and effectively.


## 5. Conclusion

This paper proposes a web app for TOF recognition. Our proposed model is termed structurally optimized stochastic pooling convolutional neural network (SOSPCNN) with explainable property achieved by Grad-CAM. The results by ten runs of 10-fold crossvalidation show that this SOSPCNN model yields a sensitivity of92.25±2.19, a specificity of 92.75±2.49, a precision of 92.79±2.29, an accuracy of 92.50±1.18, an F1 score of 92.48±1.17, an MCC of 85.06±2.38, an FMI of 92.50±1.17, and an AUC of 0.9587. Further, we develop both desktop and web apps to realize this SOSPCNN model.The shortcomings of our method are as follows: (i) our model is trained on a small dataset; (ii) our model does not go through strict medical verification; (iii) our model only considers TOF and HC.Therefore, we shall attempt to solve the above three weak points in the future. We shall try to collect more TOF and HC cardiovascular CT images. We shall invite clinicians to use our web app and return feedbacks so that we can continue to improve our model. We shall try to collect data of other heart diseases, so make our model can identify more types of diseases.



---

*Source: 5792975-2021-07-02.xml*