Keywords:

Automated Detection and Segmentation of Mitochondrial Images based on Gradient Enhancement and Adaptive Gabor Filter

Nhan Nguyen-Thanh, Tuan D. Pham and Kazuhisa Ichikawa Abstract-Information of cellular organelles location and morphology is essential for cancer simulation. In order to obtain such information, the segmentation of the organelles from electronic microscopy intracellular image is crucial. In this paper, we focus on the automatic segmentation of mitochondria organelle which is one of the most important organelles tightly related to the form of cancer. A simple three-stage strategy which includes coarse segmentation, detection and fine segmentation is proposed for fully automatic mitochondria segmentation. The local gradient calculation provides a weight factor matrix and a orientation matrix. The weight factor matrix will improve the contrast of organelle boundary of intracellular images and hence facilitate both coarse and fine mitochondria segmentation. The orientation matrix will be used for enhancing the Gabor feature extraction which make the mitochodrial detection process more accurate. Machine learning-based classifiers including knearest neighbor (k-NN), support vector machine (SVM)-based and neural network (NN)-based classifiers, are considered to learn with eight extracted features for mitochondrial detection.

Experimental results on focused ion beam (FIB) and scanning electron microscope (SEM) images of cancer cellular of human head and neck squamous cell carcinoma (SCC-61) have shown the effectiveness of proposed method.

Index Terms-cell organelles, image segmentation, image processing, mitochondria, gradient boundary enhancement, orientation adaptive, Gabor filter, SCC-61, cancer cell.

I. INTRODUCTION

In cell biology, mitochondrion is one of the main organelles that can be found in most eukaryotic cells. This type of organelle concerns with the roles of energy supply, signaling and control of cellular differentiation, growth and death. It is reported that there is a strong link between the function of mitochondria and cancer. [START_REF] Kroemer | Mitochondria in cancer[END_REF]. For example, the function of mitochondria of cancer cell is believed to be altered for resistance apoptosis [START_REF] Gogvadze | Mitochondria in cancer cells: what is so special about them?[END_REF], [START_REF] Wallace | Mitochondria and cancer[END_REF]. This leads the research to target on mitochondria for cancer therapy by changing mitochondrial metabolism or stimulating mitochondrial membrane permeabilization [4]. In addition, it is convinced that there are a relation of mitochondrial morphology and molecular controlling mechanisms during cell death [START_REF] Wasilewski | The changing shape of mitochondrial apoptosis[END_REF], [START_REF] Melino | Dynamics of mitochondrial morphology in healthy cells and during apoptosis[END_REF] and a connection of mitochondrial dynamics to key signaling cascades [START_REF] Campello | Mitochondrial shape changes: orchestrating cell pathophysiology[END_REF]. Therefore, obtaining information of location and morphology of mitochondria in intracellular space is very useful for cancer simulation.

Recent progresses on super-resolution electronic microscopy enable us to obtain this information through intracellular image. Some earlier methods such as ion-abrasion (IA) or serial block-face (SFB) can be combined with scanning electron microscopy (SEM), which is well-known for capturing surface characteristics and morphology, to obtain a lower resolution of both intercellular and intracellular space. Lately, the detail complex structure of the cell can be acquired with the combination of SEM and focused ion beam (FIB) technique, which increase the resolution of the scan. As a result, the integrated imaging system FIB-SEM has now been a powerful scientific instrument for the study of biological specimens and soft materials [START_REF] Stokes | A new approach to studying biological and soft materials using focused ion beam scanning electron microscopy (FIB SEM)[END_REF].

Automatically analyzing and understanding such kind of intracellular electron microscopy images are not a trivial task since the space of cellular contains a complicated system with variety of organelles possessing various shapes and sizes. There are some recent works that succeeds in segmentation of mitochondria [START_REF] Lucchi | Supervoxel-Based Segmentation of Mitochondria in EM Image Stacks With Learned Shape Features[END_REF]- [START_REF] Vitaladevuni | Automatic mitochondria detection in electron micrographs[END_REF] as well as some other cell organelles [START_REF] Plissiti | Combining shape, texture and intensity features for cell nuclei extraction in Pap smear images[END_REF]- [START_REF] Andres | 3D segmentation of SBF-SEM images of neuropil by a graphical model over supervoxel boundaries[END_REF]. In detail, in [START_REF] Vitaladevuni | Automatic mitochondria detection in electron micrographs[END_REF], a textural-based Gentle-Boot classifier is trained for detecting mitochondria in electron microscopy images of brain tissue. Similarly, in [START_REF] Narasimha | Automatic joint classification and segmentation of whole cell 3D images[END_REF], multiple classifiers including k-nearest neighbor algorithm (k-NN), support vector machine (SVM) and adaptive boosting algorithm (Adaboost) are used for detecting mitochondria and other organelles of highly pigmented human melanoma (MNT-1) cell based on textural feature. These works were said to be under-utilization of image information, i.e., shape information, and were improved in [START_REF] Lucchi | Supervoxel-Based Segmentation of Mitochondria in EM Image Stacks With Learned Shape Features[END_REF] by considering a graph partitioning scheme incorporated with a learn shape feature method. Though the results of these research are reasonable, the huge load works for training the system could be the most difficulty for implementation.

FIB-SEM intracellular images reveal a great detail structure of the cellular space which includes variety types of organelles such as nucleus, mitochondria, golgi apparatus, endoplasmic recticulum, cytoskeleton, vesicle, etc., with difference sizes, shapes, texture and brightness. The solid distribution of these organelles leads the segmentation task for intracellular images to be very complicated. Furthermore, in cases of images of abnormal cellulars such as cancer cells derived from malignant tumor, the structure and other characteristics of organelles tend to be more sophisticated which can cause the hard segmentation task to be harder. In this paper, images of a cancer cellular of human head and neck squamous cell carcinoma (SCC-61) [START_REF] Clark | Cortactin is an essential regulator of matrix metalloproteinase secretion and extracellular matrix degradation in invadopodia[END_REF], [START_REF] Hoshino | Establishment and validation of computational model for MT1-MMP ependent ECM degradation and intervention strategies[END_REF] captured by FIB-SEM system are investigated for segmentation mitochondria for a further cancer simulation. In order to perform this purpose, we consider a simple three-stage strategy, as shown in Fig. 1, for automatically identifying and segmenting mitochondria in intracellular images. First, an initial identification for mitochondria is implemented by coarse segmentation stage. The purpose of this coarse segmentation stage is to identify potential regions which mostly contain mitochondria. This initial identification is mainly based on the characteristic of mitochondrial shape. In other words, the output of this coarse segmentation stage are isolated areas of which the shapes are similar to the shapes of mitochondria,i.e., circular and elliptical forms derived from the cross sections of tubular mitochrodrial structures. As a result, the regions of very small and isolated organelles, which could belong to cytoplasm, free protein or small vesicle, and the regions of complex and sticky filament structures, which can originate from endoplasmic reticulum, golgi apparatus or cytoskeleton organelles, should be excluded principally. Since it is just a coarsely initial segmentation, the requirement of accuracy is not too strict. A simple combination of morphological operators are adopted to perform this task. For facilitating this coarse mitochondrial segmentation stage, a gradient boundary enhancement process, which is implemented based on a weight factor matrix formulated by the gradient magnitude of the mean transformed intracellular image, is additionally considered.

Second, the coarse segmented or candidate regions are further identified to mitochondrial and non-mitochondrial regions based on detection stage. The Gabor features and some basis features such as intensity, perimeter, and areas are extracted for mitochondrial detection stage. The gradient orientation of the mean transformed intracellular image is utilized for extracting orientation adaptive Gabor features which makes the mitochodrial detection process more accurate. In order to combine the eight proposed features for mitochondrial detection, we consider multiple types of machine learningbased classifiers including k-nearest neighbor (k-NN), support vector machine (SVM)-based and neural network (NN)-based classifiers. Since the accuracy of the three-stage segmentation method is decided in this stage, the main contribution of this

¦! "#$%&' ()$!%0!01%#)! 23#4! 56$!% "!01%#)! ())!% "!01%#)! 789@!@ &)$8 A%&B$#! Fig. 2. Mitochondrion picture [38]
paper is concentrated on this part.

Third, the detected coarse mitochondrial regions are finely segmented by fine segmentation stage. The task of this fine segmentation stage is to improve the mitochondrial segmentation by removing non-mitochondrial pixels which touch in coarse mitochondrial regions. Since the input coarse segmented regions of this stage are quite isolated and small, various kind of segmentation methods can be adopted for this stage. Herein, we adopt similar gradient boundary enhancement and combination of morphological operators, used in previous coarse segmentation stage, for performing this improving segmentation task.

The remainder of this paper is organized as follows. In section II, the coarse segmentation stage is described. In section III, the mitochondrial detection strategy is presented in detail. In section IV, the fine segmentation stage is explained. In section V, experiment results are discussed. Finally, in section VI, conclusions are drawn.

II. COARSE AND FINE SEGMENTATION

Fig. 2 illustrates structure of a mitochondrion. It can be seen that the thin outer membrane and cristae, formed by inner membranes will create thin boundaries and ridge structures in mitochondrial scans. An example of cropped FIB-SEM image of intracellular space, as shown in Fig. 4a, reveals that the space is filled by various types organelles with different sizes, shapes and brightness. These organelles are mixed together and their boundaries may not be clear. As a result, the boundary separation or boundary enhancement should be considered as a pre-processing step to facilitate the segmentation task for these organelles. Herein, we focus on thin mitochondrial boundary enhancement by adopting gradient transform.

A. Gradient transform and gradient boundary enhancement

An image I is firstly transformed such that the boundary of interested region is still preserved while the texture inside that region is suppressed. The mean filter, which is computed by

µ w (I xy ) = 1 w 2 w w I xy , (1)
where w is the window size, satisfies this requirement since it smooths texture regions while provides an average intensity

Consequently, the gradient orientation and magnitude at each pixel are computed by

ϕ xy = tan -1 ∂I xy ∂y ∂I xy ∂x (3)
and

|∇I xy | = ∂I xy ∂x 2 + ∂I xy ∂y 2 , (4)
respectively. The low intensity variation of the smoothed region causes a low gradient value while strong variation between intersection boundary and both neighbor textural region provides two new high gradient value boundaries which limit the original intersection boundary. Therefore, the value of w will determine the thickness of detected boundary. Fig. 3 illustrated the examples of gradient boundary and ridge detection. It can be seen that the gradient transform magnitude provides a good localization for boundaries between textural regions.

In order to perform boundary enhancement, gradient magnitude vector is adopted for formulating a weight factor matrix as follows:

W xy =   1 - |∇I xy | max xy (|∇I xy |)   n (5)
where n is the power transform order, which amplifies the contrast of gradient image. A larger value of n provide a higher contrast of weight factor matrix. The intensity of the original image is directly adjusted by multiplying with the weight factor matrix. In general, the value of n determine the separation level of boundary. A higher value of n can lead to more broken boundary. Therefore, the selection of n is a trade off between the separation of inter-organelle touched objects and the integrity of intra-organelle parts.

B. Mitochondrial coarse segmentation

After adopting enhancement, intracellular images is processed to obtain coarse segmented images by the following steps:

• Convert to binary images. A threshold based on Otsu method [START_REF] Otsu | A threshold selection method from gray-level histograms[END_REF] is used for getting binary version of each enhanced intracellular image. • Remove small regions. Isolated small regions in binary images which could belong to cytoplasm, free protein or small vesicle are removed.

• Remove nucleus and cell membrane region.

• Apply combination of morphological operators. The combination of morphology operators aiming to keep circular and elliptical shape is described as follows:

-Fill holes, -Apply opening operator with disk structuring element, -Remove small regions, -Compute erosion by disk structuring element, -Remove small regions, -And compute dilation by disk structuring element. Since the morphological operators are adopted to a large area with various condition, the morphological operators is principally selected such that the integrity of mitochondrial organelles is warranty. This selection means that a level of merged objects is acceptable.

C. Mitochondrial fine segmentation

The fine segmentation stage is executed for finely segmenting detected coarse segmented mitochondrial regions. Due to the separation of detected coarse mitochondrial regions, various segmentation method can be applied to provide a fine segmented output without affecting other areas. Therefore, the fine segmentation stage is totally opened for multiple segmentation methods.

Herein, a similar gradient boundary enhancement and morphological operator combination, used in previous coarse segmentation stage, are performed for finely segmenting mitochondria. However, unlike coarse segmentation stage, of gradient boundary enhancement and morphological operators are only affect in small detected mitochondrial regions. Therefore, the parameters are selected such that a good fine segmentation result can be obtained.

Two examples of fine segmentation results are illustrated in Fig. 5c and5f. It can be seen that the fine segmentation results are largely improved compared to coarse segmented versions of these objects which are shown in Fig. 5a and5b.

III. DETECTION

A. Feature extraction

Extracting an efficient feature set is a crucial for any detection algorithm. A good feature set will provide a good classification result. In this work, the feature set is utilized for determining whether a coarse segmented region contains mitochondria. The best feature for discriminating a mitochrondion from other types of organelles is the ridge texture feature formed by the cross section of cristae. In addition, since the coarse segmented regions can be any types of organelles, some basic features such as intensities, perimeter, area, etc, can provide a certain accuracy level for detecting mitochondria. As a result, we select the feature set that specifies a mitochondrial region mainly based on texture features and some basic regional properties as shown in Fig. 6.

1) Orientation adaptive Gabor filter and Gabor features:

Gabor filter is selected for extracting textural feature which objects to ridge structure of mitochondria. The reason for choosing Gabor filter is due to the fact that Gabor filtering is well-known as an invariant and effective method for extracting textural features [START_REF] Grigorescu | Comparison of texture features based on Gabor filters[END_REF]- [START_REF] Bianconi | Evaluation of the effects of Gabor filter parameters on texture classification[END_REF]. Furthermore, Gabor filter is convinced to be well-adapted to fiber or ridge textural extraction which was contained in variety applications such as fingerprint classification [START_REF] Zhang | Spectral analysis based fingerprint image enhancement algorithm[END_REF], [START_REF] Gonzaga De Oliveira | A methodology for identification of fingerprints based on Gabor filter[END_REF], Iris recognition [START_REF] Chung-Chih | Iris Recognition Using Possibilistic Fuzzy Matching on Local Features[END_REF]- [START_REF] Sanchez-Avila | Two different approaches for iris recognition using Gabor filters and multiscale zero-crossing representation[END_REF], prostate cancer image segmentation and recognition [START_REF] Dinggang | Segmentation of prostate boundaries from ultrasound images using statistical shape model[END_REF], [START_REF] Mohamed | Prostate Cancer Spectral Multifeature Analysis Using TRUS Images[END_REF], blood vessel detection [START_REF] Wu | On the adaptive detection of blood vessels in retinal images[END_REF], character recognition [START_REF] Wang | Gabor filters-based feature extraction for character recognition[END_REF], [START_REF] Jonghyon | Online signature verification using temporal shift estimated by the phase of Gabor filter[END_REF]. This effectiveness can be explained by the property of both frequency and orientation selective which provides a joint optimal resolution in both frequency and spatial domain. The general 2D Gabor filter function is defined by

g λ,θ (x, y) = exp - 1 2 x 2 θ σ 2 x + y 2 θ σ 2 y exp j 2π λ x θ (6)
where x θ = x cos θ + y sin θ and y θ = -x sin θ + y cos θ, λ is the wavelength of the sinusoidal plane, θ is the orientation of the Gabor filter, and σ x and σ y are the standard deviations of the Gaussian envelope along the x axis and y axis, respectively. The Gabor filter can be analyzed into even and odd parts by g λ,θ = g even λ,θ + jg odd λ,θ . Denote γ = σ y /σ x and σ = σ x = γσ y then the odd and even parts of the Gabor filter are g odd λ,θ,γ (x, y) = exp -

x 2 θ + γy 2 θ 2σ 2 cos 2π λ x θ (7)
and

g even λ,θ,γ (x, y) = exp - x 2 θ + γy 2 θ 2σ 2 sin 2π λ x θ (8)
respectively. The half-response spatial frequency bandwidth B is selected by [START_REF] Grigorescu | Comparison of texture features based on Gabor filters[END_REF]:

B = log 2 σ λ π + ln 2 2 σ λ π -ln 2 2 (9)
The Gabor filter after defined will be convoluted with the input image I(X, Y ). The magnitude of Gabor coefficient at a pixel is computed as follows:

G θ (X, Y ) = w 2 x=-w 2 w 2 y=-w 2 I (X + x, Y + y) g λ,θ,γ (x, y) (10
) where I is a window of size w × w.

In general, the Gabor coefficient of a pixel can be calculated as a square root combination of coefficients determined by multiple Gabor filters with different orientations. However this computation method will decrease the regional contrast of output coefficients since "out-of-orientation" Gabor coefficients will cause a reduction on the "in-orientation" ones which is usually the biggest Gabor coefficient. Furthermore, high-oriented ridge texture is frequently appeared in mitochondrial regions, while do not exist in other organelles' regions. Therefore, determinating a correct orientation θ for Gabor filter is necessary for enhancing the discrimination between mitochondrial and non-mitochondrial regions.

In this paper, we propose an orientation adaptive Gabor filter for such textural feature extraction. Output Gabor coefficient of a pixel is computed by a Gabor filter of which the orientation θ xy is calculated based on to the pixel's local gradient orientation ϕ xy in (3) as follows:

θ xy = ϕ xy + π 2 . ( 11
)
In order to reduce the computation of the orientation adaptive Gabor filter, the orientation for Gabor filter at each pixel θ xy is linearly quantized to a certain set of N orientations {θ k } N k=1 as follows:

θ ′ xy = kπ N (k-1)π N < θ xy ≤ kπ N , k = 1, ..N -1 0 θ xy > (N -1)π N ( 12
)
Due to the the matching of orientation between Gabor filter and mitochodrial ridge texture, the maximum value of Gabor coefficients of a mitochodrial region will be higher than that of other non-mitochondrial regions. Furthermore, since most of textures of non-mitochondrial regions have more unvarying distribution than that of mitochodrial regions, the standard deviation of Gabor coefficients of mitochondria regions tends to be higher than that of non-mitochondrial ones. For that reason, we consider the maximum and the standard deviation values of Gabor coefficients are the most distinctive features for detecting mitochondria.

In order to implement Gabor features extraction for mitochondiral detection, as shown in Fig. 6, we adopt adaptive Gabor filters to the high pass filtered version and the gradient boundary enhanced version of intracellular images. The reason for this utilization can be explained by the fact that the textural regions of mitochondria in high pass filtered version and gradient boundary enhanced version are more highlighted than that in the original intracellular image.

2) Basis regional features: Complementing to textural features, we utilize some other basic quantitative parameters for the feature set. For distinguishing mitochondrial region from some darker or brighter regions of other organelles such as vesicles, cytoskeleton, etc, the mean value of the intensity of coarse segmented region is added to the feature set.

In addition, since ridges and boundaries of mitochondria, formed by brighter pixels compared to the remainders, are thinner than those of other organelles, the ratio of mean original intensity of pixels belonging to Otsu thresholded gradient boundary enhanced coarse segmented regions divided by that belonging to Otsu thresholded high pass filter coarse segmented regions is extracted as a feature for distinguishing mitochondrial and non-mitochondrial regions. We denote this term as M eanHP F/M eanGrad. in Fig. 6.

Another feature of the regional ratio of perimeter divided by area, which can help to recognize the filament structures derived from endoplasmic reticulum, golgi apparatus or cytoskeleton organelles, is also utilized. Normally, the ratio of regional perimeter and area of a mitochondrion is higher than that of a filament region. We signify this ratio as P erimeter/Area in Fig. 6.

Besides, it is obvious that the solidity, which is defined by ratio of binary regional area divided by convex area, of a mitochondrial region is differ from the solidity of most of other organelles' regions, except some special touched filament regions. Furthermore, since the filled area of a mitochondrial region is almost equivalent to its convex area while the filled areas of some of special touched filament regions differ from their convex areas, the replacement of convex area by filled area will increase the distinction of the solidity feature between mitochondrial regions and some special touched filament regions. Therefore, the modified solidity which is re-defined by the ratio of binary regional area divided by filled region area, is also utilized for detecting mitochondrial regions. We indicate this modified solidity as Area/F illedArea in Fig 6.

B. Machine learning-based mitochondrial detection

In order to detect mitochondrial regions, we consider multiple types of classifier including k-NN, SVM and neural network (NN). The coarse segmented regions is divided into a training set and a testing set. Extracted feature set of the training set is combined with validation data in a learning process. The output system of that learning process, which will be used for detecting mitochondrial regions, is evaluated based on the testing set.

It should be noticed that the proposed feature set for mitochondrial detection is only contained in an 8-dimensional space. Therefore, the resources for training task of the proposed method is much lower than those of many of previous works [START_REF] Lucchi | Supervoxel-Based Segmentation of Mitochondria in EM Image Stacks With Learned Shape Features[END_REF], [START_REF] Narasimha | Automatic joint classification and segmentation of whole cell 3D images[END_REF]. However, the proposed feature set still ensure good coverage for detecting mitochondrial regions, since both basis geometrical features and textural features of coarse segmented regions are utilized. Therefore, the proposed method is simple in machine learning implementation aspect.

IV. EXPERIMENTS

The experiments are conducted with FIB-SEM intracellular images collected from the cancer cell line of a human head and neck squamous cell carcinoma (SCC-61). Fig. 4a shows an example FIB-SEM scan of a cancer cell SCC-61. Compared to the images of CA1 hippocampus and striatum tissues [START_REF] Lucchi | Supervoxel-Based Segmentation of Mitochondria in EM Image Stacks With Learned Shape Features[END_REF] which are also captured by FIB-SEM system, the studied images images are more complicated. Moreover, due to higher resolution, the mitochodrial segmentation task for the studied FIB-SEM images of SCC-61 cancer cell are more difficulty than the mitochodrial segmentation task for IA-SEM images of MNT-1 cells [START_REF] Narasimha | Automatic joint classification and segmentation of whole cell 3D images[END_REF], [START_REF] Heymann | 3D Imaging of mammalian cells with ion-abrasion scanning electron microscopy[END_REF], [START_REF] Murphy | Correlative 3D imaging of whole mammalian cells with light and electron microscopy[END_REF] and SBF-SEM images of mouse cerebellum, dentate gyrus, and CA3 hippocampus [START_REF] Giuly | Method: automatic segmentation of mitochondria utilizing patch classification, contour pair classification, and automatically seeded level sets[END_REF].

In general, the performance of the proposed three-stage mitochondrial segmentation method is mainly decided in the accuracy of the detection stage. Therefore, we focus on analyzing the effectiveness of the proposed detection method in this section.

A. Parameters selection for feature extraction

There is a resonance effect of spacial frequencies of mitochondrial ridge texture and Gabor filter. As a result, the selection of spacial frequency for Gabor filter is very important for making an accurate mitochondrial detector. Herein, we adopt Fourier transform-based frequency analysis method for determining spacial frequency λ of ridge textures of mitochondria. The determination process is illustrated in Fig. 7. A slice cut of a single mitochondrial region, e.g., the solid line as shown in Fig. 7a, is extracted. Fig. 7b plots the spacial intensity series of the slice cut in Fig. 7a. Adopting Fourier transform to the spacial series, we obtain the magnitude spectrum of the spacial series as shown in Fig. 7c. The spacial frequency of the ridge texture is determined to be equal to the first order harmonic component, i.e., the peak at λ = 5 in Fig. 7c.

It should be noted that the slice cut n Fig. 7a is intentionally selected such that the cut is orthogonal with the ridge texture orientation. In case of a cut slice is "out-of-orientation" with ridge texture, the output frequency may be distorted. Consequently, the output of Gabor filter will strongly depend on the local orientation of mitochondrial ridge texture. This effect also shows the validity and potential of the proposed idea on using Gabor filter with adaptive orientation.

Obviously, the proposed adaptive Gabor filter is designed for resonant with mitochondrial ridge texture. However, a coarse mitochondrial region can also contain some touched parts of other organelles, which can create low values of Gabor coefficient. As a result, the best Gabor feature for mitochondrial detection is the maximum value of Gabor coefficients of coarse segmented region. Next, we investigate the influence of parameters on the mitochondrial detection efficiency of the key detection factor of maximum value of regional Gabor coefficients. Hereafter, we consider high pass filtered version of original intracellular space scans for the investigation. The results are shown in Fig. 8. In order to present the effectiveness of detection parameters, we utilize receiver operating characteristic (ROC) curves, which is the plot of true positive rate (TPR) and false positive rate (FPR). In this paper, TPR is defined by the ratio of number of correct detected mitochondrial regions divided by the total number of mitochondrial regions in the images. And FPR is defined by the ratio of number of incorrect detected mitochondrial regions divided by the total number of non-mitochondrial regions.

Fig. 8a presents ROC curves of the maximum value of coarse segmented regional Gabor coefficients computed based on the proposed orientation adaptive Gabor filter corresponding to different spacial frequencies when the parameter of window size w and the orientation quantization levels N are fixed at 5 and 24, respectively. As shown in Fig. 8a, the mitochondrial detection ability of maximum regional Gabor coefficient computed by the proposed method depends on the value of spacial frequency. Matching with the previous analysis, the detection performance of maximum regional Gabor coefficient at spacial frequency of λ = 5, which is equal to the first order harmonic frequency component, is higher than that at other frequencies. This observation proves the effectiveness of choosing correct spacial frequency for mitochondrial detection.

Since for our proposed orientation adaptive Gabor filter, the calculation of local orientation of intracellular images strongly affect to the output Gabor coefficients and hence the detection performance of the maximum Gabor coefficient feature. There are two parameters that influence on the calculation of local orientation of intracellular images. They are the window size w for computing the local gradient orientation and the number of orientation quantization levels.

First, the impact of window size w to the detection performance of the maximum Gabor coefficient feature is presented on Fig. 8b. From the ROC curves corresponding to three value of w of 5, 7, and 9, it can be seen that the value of w do not highly affect on the mitochondral detection performance of the maximum Gabor coefficient feature since the ROC curves are almost the same.

Second, we investigate the effect of number of orientation N on mitochondrial detection performance of the maximum Gabor coefficient features calculated by both the proposed N quantized orientation adaptive and the conventional root mean square N orientations Gabor filter. The results are plotted in Fig. 8c. It can be seen that the accuracy of mitochondrial de-tection of the maximum Gabor coefficient features of the both Gabor filter will slightly increase along with the increment of orientation number N .

Particularly, it can be also seen in Fig. 8c that the mitochondrial detection accuracy of our proposed orientation adaptive Gabor filter is higher than that of conventional Gabor filter in all case of N . The proposed oriental adaption for Gabor filter following local orientation of mitochondrial ridge texture improves the detection ability for Gabor filter.

Finally, in order to confirm the effectiveness of Gabor filter for mitochondrial detection, we compare the mitochondrial detection ability of our proposed Gabor filter-based textural features including the maximum and standard deviation values of coarse segmented regions with other textural features including homogeneity, contrast, correlation and entropy derived from famous gray level co-occurrence matrix (GLCM). The results is presented in Fig. 8d. Obviously, our proposed features of maximum and standard deviation values of coarse segmented regions outperforms the textural features extracted based on GLCM for mitochondria detection.

B. Machine learning-based mitochondrial detection

Our experiments was conducted on 20 scans of a SCC-61 cell. Half of the images are used for training and the remainders are used for evaluating the trained models. Multiple machine learning methods including k-NN-based, SVM-based and NN-based classifier for mitochondria detection are considered. For k-NN classifiers, we tested with multiple values of k including k = 5, 9 and 15. For SVM-based classifiers, sequential minimal optimization (SMO) [START_REF] Platt | Sequential Minimal Optimization: A Fast Algorithm for Training Support Vector Machines[END_REF] and least square (LS) methods are considered for training SVMs. For NN-based classifiers, we performed the training with four configurations of feed-forward back-propagation neural networks, i.e., two configurations with one hidden layer including 10 nodes (net-1) and 20 nodes (net-2) and two configurations with two hidden layers including 10 and 3 nodes (net-3), and 20 and 8 nodes (net-4).

The values of the spacial frequency λ and the number of orientations N of Gabor filters for feature extraction are selected as 5 and 16, respectively, which are reasonable values as analyzed in previous subsection. The value of window size of w, which has shown a low influence to the detection Next, we keep w = 7, n = 2, and λ = 5 and change the value of N in order to investigate the affect of orientation number to final detection performance. As shown in Table II, there is a saturation of final detection accuracy (ACC) when N larger than 8. The values of TPR and FPR in such cases are maintained as about 92% and 16%, respectively. These value means that the error rate in identifying mitochondria regions and non-mitochondrial regions are corresponding to about 8% and 16%.

Most of the cases that cause these inaccurate identification belong to the small mitochondrial regions which are falsely merged with large neighbor organelles in coarse segmentation stage, and small non-mitochondrial regions of which the texture can not be detectable. This observation is demonstrated in Fig. 9 where the Fig. 9a indicates coarse segmentation regions of an intracellular image, and Fig. 9a expresses detected coarse mitochondrial regions obtained by an SVM-based mitochondrial detector with λ = 5, γ = 0.5, B = 1, w = 7, n = 2 and N = 16. In detail, the small mitochondrial regions that are missed due to the affect of large merged neighbor organelles are annotated by (A) and (B). The small non-mitochondrial region that are falsely detected because of the similarity of its texture to mitochondrial small region's texture is marked by (C).

Since most of incorrect detected and miss detected regions are small, the accuracy of our proposed algorithm computed in pixel will be higher than 90%. Therefore, the accuracy of our method can be comparable with the results in [START_REF] Lucchi | Supervoxel-Based Segmentation of Mitochondria in EM Image Stacks With Learned Shape Features[END_REF], [START_REF] Narasimha | Automatic joint classification and segmentation of whole cell 3D images[END_REF], i.e., about 95%, while our proposed algorithm is much simpler. Furthermore, it is possible to recovered or discarded correctly the small incorrect detected and miss detected regions by applying 3D reconstruct techniques to adjacent slices where these regions appear clearer to be classified correctly.

V. CONCLUSIONS

We have presented a simple three-stage method for mitochondrial segmentation. The proposed method initially identify mitochondria by coarse segmentation and further identify by a machine learning-based detector. A gradient boundary enhancement is utilize for facilitating the coarse segmentation. A orientation adaptive Gabor filter is adopted for effectively extracting features for mitochondrial detection. The proposed method has illustrated an appropriately simple and high accuracy implementation for the hard problem of segmentation mitochrondria in intracellular space.

Fig. 1 .

 1 Fig. 1. Mitochondrial segmentation strategy

Fig. 4

 4 shows a comparison example of coarse mitochondrial segmentation outputs for different enhancement types. Fig 4b is the coarse segmentation output of the original intracellular image in Fig 4a. Fig 4c exhibits the enhanced version based on contrast limited adaptive histograph equalization (CLAHE) [37] of the original intracellular image. Fig 4d gives the coarse segmentation regions obtaining by adopting the proposed algorithm to Fig 4c. Similarly, Fig 4f displays the coarse segmentation regions of the gradient boundary enhanced intracellular image shown in Fig 4e.It can be seen that for the same coarse segmentation algorithm, the gradient boundary enhanced image provide a better coarse segmentation result compared to original non-enhanced image and CLAHEbased enhanced image. Furthermore, since mitochondria are contained in most of coarse segmented regions, a further identification for mitochondrial regions will be implemented in the next detection stage.

Fig. 4 .

 4 Fig. 4. Coarse segmentation comparison examples. (a) An original intracellular image. (b) Coarse segmentation of (a). (c) Enhanced image of (a) based on CLAHE [37]. (d) Coarse segmentation of (c). (e) Enhanced image of (a) based on gradient boundary enhancement. (f) Coarse segmented image of (e).

Fig. 5 .

 5 Fig. 5. Fine segmentation examples. (a) Detected coarse mitochondrial segmentation object 1. (b) Retrieved original cropped region containing object 1. (c) Fine segmented of object 1. (d) Detected coarse mitochondrial segmentation object 2. (f) Retrieved original cropped region containing object 2. (g) Fine segmented of object 2.

Fig. 6 .

 6 Fig. 6. Feature extraction scheme.

Fig. 7 .

 7 Fig. 7. Spacial frequency determination example. (a) A mitochondrial region with a horizontal slice cut. (b) Intensity of pixels on the slice cut. Fourier transform of intensity series(b).

Fig. 8 .

 8 Fig. 8. ROC curves of regional maximum Gabor coefficient in comparison of (a) multiples space frequency λ values, (b) conventional Gabor filter and orientation adaptive Gabor filter, (c) different values of window size w for determining orientation, (d) of other Gabor and gray level co-occurence matrix (GLCM) features.

TABLE II AFFECT

 II OF ORIENTATION NUMBER ON DETECTION ACCURACY.

		ACC	TPR	FPR
	N = 4	87.32	91.91	20.56
	N = 8	89.37	92.88	16.67
	N = 16	90.18	92.88	14.44
	N = 24	88.55	92.23	17.78
	N = 32	89.16	92.56	16.67

accurate classifier with approximate 90% of correct detection.

ACKNOWLEDGMENTS

We thank Dr. Daisuke Hoshino and Professor Motoharu Seiki of Division of Cancer Cell Research, Institute of Medical

performance of the key feature of maximum regional Gabor coefficient extracted from high pass filtered versions of intracellular scans, are set with four different values of 3, 5, 7, and 9 to consider its effect to the final detection accuracy. Similar to w, the values of power order are also selected as 2 and 3 to evaluate its impact to the extracted features of maximum and standard deviation values of regional Gabor coefficients derived from gradient boundary enhanced versions of intracellular scans.

The detection accuracy, which includes true positive and true negative rates, of multiple investigated classifiers are shown in Table I. It can be seen that the influence on final detection performance of window size w value is higher than that of power order n value. This observation can be explained by the fact that the value of w affect to both high pass filtered and gradient boundary enhanced versions of original scans while the value of n only have an effect on the gradient boundary enhanced ones. Similar to spacial frequency λ, there is also a better choice of window size w where the detection performance is higher, i.e., w = 7, in this experiment.

For all of investigated cases presented in Fig. I, it can be seen that at least 76% of coarse segmentation regions are classified correctly. Obviously, the method of classifier do influence on the accuracy of detection, i.e., the detection performance generally decreases in the order of SVM-based, NNbased and k-NN classifier. For training selection, it is favorable to choose SVM-based classifier trained by LS method since the detection accuracy of this SVM method is almost the highest Science, University of Tokyo, for providing cancer cell images.


Nicolas Martin

Nicolas Mathieu

Nico Pallamin

Martin Ragot

Jean- Marc Diverrez

AUTOMATIC RECOGNITION OF VIRTUAL REALITY SICKNESS BASED ON PHYSIOLOGICAL SIGNALS

Virtual Reality (VR) sickness seems one of the main limitations to the large-scale adoption of VR technologies. This disturbance seems to induce physiological changes that affect the sympathetic and parasympathetic activities of the users. Thereby, it seems relevant to measure users' physiological data in order to prevent and reduce VR sickness. This paper presents the results of an initial real-life experiment of VR sickness detection based on physiological data. The electrodermal, cardiac and subjective data of 27 participants was recorded during VR sessions. Machine Learning algorithms were trained and the best model (Gradient Boosting) explained 48% of the VR sickness variance. These results demonstrate the opportunity to develop an automatic and continuous tool to detect the appearance of VR sickness based on physiological signals. This tool will prove very valuable to the VR industry.

INTRODUCTION

Virtual Reality (VR) appears as a major technological breakthrough and a main business opportunity for the entertainment industry. The VR market is expected to expand exponentially with worldwide revenues for the AR/VR growing to more than $162 billion in 2020 [START_REF] Idc | Worldwide Semiannual Augmented and Virtual Reality Spending Guide[END_REF]. However, one main limitation to its large-scale adoption is VR sickness especially because of health, ethical, legal and acceptability aspects. VR sickness is a common problem that could affect up to 60% of adult users [START_REF] Chang | Postural activity and motion sickness during video game play in children and adults[END_REF].

The necessity to better detect and prevent the appearance of VR sickness is at the origin of this research cooperation between b<>com and Editorial user research lab, Ubisoft Paris aiming at the development of an automatic VR sickness detection tool.

VIRTUAL REALITY SICKNESS: CAUSES AND SYMPTOMS

Previous research proposed several hypotheses to explain VR Sickness. The most common theory is the sensory conflict theory [START_REF] Reason | Motion sickness adaptation: a neural mismatch model[END_REF]. According to this theory, VR sickness is the result of conflicts and inconsistency between the different sensorial information sent to the brain when the user evolves in the virtual world. In VR, the most common conflict is a discrepancy between the motion information coming from the vestibular system and the motion information coming from the visual system. Indeed, this information comes from two separate systems: the vestibular system located in the inner ear and the visual system [START_REF] Mousavi | A Review on Cybersickness and Usability in Virtual Environments[END_REF]. This discrepancy will be detected by the brain and will induce, in sensitive participants, the symptoms of VR sickness.

VR sickness is a complex phenomenon and while motion cues play a primary role, multiple factors are known to contribute to the appearance of the sickness. Three main categories of factors can be identified:

1) The factors related to the characteristics of the stimuli: the spatial frequency [START_REF] Diels | Frequency characteristics of visually induced motion sickness[END_REF], the reactivity of the system [START_REF] St | The effects of 0.2Hz varying latency with 20-100ms varying amplitude on simulator sickness in a helmet mounted display[END_REF], the wideness of the Field Of View (FOV) [START_REF] Fernandes | Combating VR sickness through subtle dynamic field-ofview modification[END_REF], etc.

2) The factors related to the predispositions of the users: gender [START_REF] Paillard | Motion sickness susceptibility in healthy subjects and vestibular patients: effects of gender, age and trait-anxiety[END_REF], age [START_REF] Arns | The relationship between age and incidence of cybersickness among immersive environment users[END_REF], the predisposition to suffer from migraine attack [START_REF] Cuomo-Granston | Migraine and motion sickness: what is the link?[END_REF], etc.

3) The users' past experiences [START_REF] Bertolini | Moving in a Moving World: A Review on Vestibular Motion Sickness[END_REF], according to the sensory conflict theory. In fact VR sickness appears only when a present set of vestibular and visual information is not congruent with what is expected from previous experiences [START_REF] Oman | A heuristic mathematical model for the dynamics of sensory conflict and motion sickness[END_REF].

The complex nature of VR sickness is not limited to its different causes, but it is reflected in its symptomatology as well. While the most evident and detrimental symptom of VR sickness is the nausea, the complete symptomatology of VR sickness includes other elements like general discomfort, headache, disorientation and eye strain. The intensity as well as the duration of the symptoms are quite variable. They depend on the characteristic of the stimulus and the user predisposition to VR sickness. In the majority of cases, the symptoms disappear some minutes after the end of the stimulation. Nevertheless, there are documented cases when the symptoms were still present 6 hours after the VR experience [START_REF] Kolasinski | Simulator sickness in virtual environments[END_REF].

THE ASSESSMENT OF VIRTUAL REALITY SICKNESS

The traditional way to evaluate VR sickness is based on subjective questionnaires. Various questionnaires have been previously developed to assess a sickness level. The most popular is probably the "Simulator Sickness Questionnaire" (SSQ) proposed by Kennedy in 1993 [START_REF] Kennedy | Simulator Sickness Questionnaire: An Enhanced Method for Quantifying Simulator Sickness[END_REF]. It is constituted of 16 questions to evaluate 3 categories of symptoms: nausea related symptoms, oculomotor symptoms and disorientation symptoms.

While still widely used for the simplicity of its deployment and analysis, the SSQ has several limitations. Indeed, this measure is global, gathered a posteriori, intrusive and punctual. These limitations make the SSQ (and any other questionnaire) inadequate to be integrated in an automatic tool to evaluate VR sickness. Another method was proposed by Keshavarz and collaborators [START_REF] Keshavarz | Validating an efficient method to quantify motion sickness[END_REF] which consists in requesting the participants to verbally give their evaluation of sickness level on a scale ranging from 0 to 20, with zero representing no discomfort and twenty representing barely supportable sickness. This approach seems to be a very good alternative to the SSQ since very strong correlation (r = .79) between verbal rating and SSQ scores were reported in this study. Nevertheless, this method is cognitively disruptive for the user experience in a context of VR-video games.

To be effective and valuable for the VR industry, the assessment method must ideally be continuous, automatic and real time. Considering these constraints, physiological data related to the perceived VR sickness was used in this paper.

The physiological measures of VR sickness

The physiological response associated to VR sickness is due to the connections between the vestibular and the autonomic nervous system. The conflicting inputs from visual, vestibular and somatosensory afferents induce a vestibular autonomic response. This response involves both the sympathetic and the parasympathetic systems and affects for instance the heart rate variability and the skin conductance [START_REF] Patrão | How to Deal with Motion Sickness in Virtual Reality[END_REF].

Various researches support the hypothesis of a measurable VR sickness by physiological response. For instance, Gavgani and collaborators [START_REF] Gianaros | Relationship between temporal changes in cardiac parasympathetic activity and motion sickness severity[END_REF] underline a correlation between phasic skin conductance activity and the reported nausea ratings. While Dennison [START_REF] Dennison | Use of physiological signals to predict cybersickness[END_REF] found relations between subscale scores of the SSQ and bradygastric power, breathing rate, pulse amplitude and blinking. Moreover, Ohyama and collaborators [START_REF] Ohyama | Autonomic responses during motion sickness induced by virtual reality[END_REF] investigated the cardiac responses associated to VR sickness. But, cardiac responses didn't statically correlate with the user subjective evaluation.

These findings suggest the existence of a relation between the physiological responses and VR sickness. Nevertheless, such relation seems too complex to be exploited using classical analyses. To handle such complex data, Machine Learning approach has spread in recent years to go beyond the limits of classical statistical approaches. Nam and collaborators [START_REF] Nam | Automatic detection of nausea using bio-signals during immersion in a virtual reality environment[END_REF] use such method to detect nausea in VR using various bio-signals (EEG, ECG, PPG, SCL, EOG, SKT). Results seem promising as they are able to partially detect nausea in real time [START_REF] Nam | Automatic detection of nausea using bio-signals during immersion in a virtual reality environment[END_REF].

Automatic recognition of VR sickness: the Machine Learning approach

To automatically detect VR sickness using physiological data, the use of Machine Learning algorithms seem relevant. In this approach, physiological and declared VR sickness are mapped using Machine Learning algorithms. In others words, the training of Machine Learning algorithms, using the supervised techniques [START_REF] Kotsiantis | Supervised Machine Learning: A Review of Classification Techniques[END_REF], aims to infer the function between the input data (i.e., physiological data) and the output data (i.e., subjective labels). Indeed, Machine Learning algorithms have the ability to learn without being explicitly programmed [START_REF] Bell | Machine Learning: Hands-On for Developers and Technical Professionals[END_REF]. After training, the models should be able to automatically recognize in real time, for any new user, the VR sickness related to the physiological data without requesting a subjective response.

Physiological signals and subjective VR sickness level were collected to investigate the possibility of detecting VR sickness automatically. A large and relevant labeled dataset (i.e., physiological data labeled with subjective evaluation of VR sickness) to train models is fundamental to obtain an efficient system of VR sickness detection. So, Editorial user research lab, Ubisoft Paris and b<>com decided to cooperate with the intent to collect and analyze consistent amounts of physiological data of test users viewing VR content. Contrary to previous studies [START_REF] Ohyama | Autonomic responses during motion sickness induced by virtual reality[END_REF], [START_REF] Nam | Automatic detection of nausea using bio-signals during immersion in a virtual reality environment[END_REF] the aim of the current paper is to develop and evaluate a new solution more adapted to the industrial context than previous approaches based on invasive and/or expensive sensors (e.g., EEG or eye tracking).

METHOD

27 participants (22 men and 5 women -average age: 25.93 years; standard deviation: 4.39) were recruited by Ubisoft. Participants were requested to test three VR game prototypes during 30 minutes. The prototypes were interactive games currently under development and they implied a relevant amount of camera's movements. The first game was a space simulation using a third-person point of view. The second was an arcade car racing game (in first-person and third-person point of view) and the third was a space firstperson shooter. To increase the probability of VR sickness we didn't implement in these prototypes any of the typical countermeasures (i.e. Reduction of the peripheral vision) applied in the gaming industry to prevent VR sickness. The participants were informed that they were free to stop the experiment at any time.

In order to induce various levels of VR sickness, the game's levels tested by the participants were designed to increase progressively the VR sickness induction. VR content was presented using the "Oculus Rift"1 ® or "HTC Vive"2 ® devices.

Measures

Two types of measures were collected: the physiological data and subjective data of the users' sickness.

The physiological data was collected using the "Shimmer GSR+"®3 . This sensor measures the electrodermal activity (EDA) and the Blood Pulse Volume (BVP) using PhotoPlethysmoGraphy (PPG) as the measure of cardiac activity. EDA sensors were placed on the middle phalanx of the non-dominant hand two first fingers. The data was collected in a continuous way during the experimental session and the sampling frequency was set at 128 Hz.

The subjective data was collected following the method proposed by Keshavarz and collaborators [START_REF] Keshavarz | Validating an efficient method to quantify motion sickness[END_REF] presented in the previous section. Participants were instructed to express their sickness level evaluation on a scale from 0 to 20 in response of auditory stimuli presented every 40 seconds. To grant the accuracy of the reply and the correct time alignment, the audio of the responses was recorded during each session and transcribed at the end of each experimental session.

Features extraction

The raw physiological signals appear very complex to interpret directly. So, it is generally necessary to extract specific features from raw signals. To extract these features, a toolbox was specially developed in Python. For this purpose, the raw signals were divided into subsamples of 30 seconds. For each of the 715 reported subjective responses, the last 30 seconds of physiological data (BVP and EDA) are used to extract the physiological features.

Concerning the cardiac activity, two types of features are extracted [START_REF] Kreibig | Autonomic nervous system activity in emotion: a review[END_REF]: time-domain features (HR, AVNN, SDNN, nn50, pNN50) and frequency-domain features (VLF, LF, LF/HF, Total Power).

Concerning the EDA, the signal is composed of two components: a phasic and a tonic part [START_REF] Boucsein | Electrodermal Activity[END_REF]. The phasic part (also called Skin Conductance Level -SCL) corresponds to slow change in the EDA while the tonic part (also called Skin Conductance Responses -SCR) corresponds to the rapid physiological responses to a stimulus. The process of extracting those two components from the raw signal is constituted of different steps. A low-pass filter (cutoff frequency = 1Hz, order = 3) is firstly applied to reduce noise in the signal. Then, the tonic part of EDA is extracted from this filtered signal using a low-pass filter (cutoff frequency = 0.05Hz, order = 3) [START_REF] Braithwaite | A guide for analysing electrodermal activity (EDA) & skin conductance responses (SCRs) for psychological experiments[END_REF] and averaged over the whole filtered signal (mean EDA). The phasic part is obtained by subtracting the tonic signal to the filtered signal. On the phasic part, the SCR are manifested by peaks on signals. On these peaks 4 , the amplitude (mean, standard deviation, minimum and maximum) and half of recovery time (mean, standard deviation, minimum and maximum) are measured. Moreover, an Empirical Mode Decomposition [START_REF] He | Multimodal Affective Dimension Prediction Using Deep Bidirectional Long Short-Term Memory Recurrent Neural Networks[END_REF] was applied to EDA signal to extract Intrinsic Mode Functions (IMF).

RESULTS

Analysis of correlation

Based on the extracted features, Pearson correlations between the physiological features and the subjective data were computed. Only five correlations appear as significant. Of these five correlations, three concern the EDA features: the Mean Amplitude of the EDA (r = .12), the Minimal Amplitude of the EDA (r = -.14) and the Maximal Amplitude of the EDA (r = -.07). The two remaining correlations concern cardiac features: the Heart Rate (r = .22) and the NN50 (r = .18).

The majority of the 33 computed correlations appear as not significant. This can be the consequence of the multidimensional and nonlinear relationships between the physiological features and the subjective evaluation of VR. So, to handle such complex data and recognize VR sickness, Machine Learning techniques were used.

Machine Learning

In order to create an automatic tool to recognize VR sickness, Machine Learning models were trained. To test the performance of trained model, the dataset was randomly split in two independent parts: a training dataset and a testing dataset corresponding respectively to the 80% (572 samples) and the 20% (143 samples) of the original dataset. Then, three Machine Learning models were trained to select the most efficient model for our task. The three selected models were: 1) Random Forest, 2) Gradient Bosting and 3) Support Vector Machine.

Regression models were chosen instead of classification models. Indeed, a continuous measure of VR sickness seems more adapted to represent the human perception compared to exclusive classes. Three metrics were computed for each method: RMSE (Root Mean Square Error), MAE (Mean Absolute Error) and R-Square (i.e., coefficient of determination) 5 . According to these results, Gradient Boosting method seems to offer the best result with more variance explained and smaller error.

CONCLUSION

VR is potentially a major business opportunity with an exponential growth. But serious acceptability issues, in particular VR sickness, could interfere with the uptake of this technological breakthrough by the mass market. In this paper, a solution to detect the development of VR sickness based on physiological signals was tested. Models were trained to recognize VR sickness using Machine Learning methods and showed promising results (the best model explained 48% of the VR sickness variance). Indeed, the ability to evaluate in a continuous and automatic way the appearance of the users' sickness seems a necessary stage to prevent VR sickness. Further research work will be conducted to confirm the results and evaluate the performance of the tool in real time.

Lastly, this type of automatic system able to assess the level of VR sickness of the end user seems extremely valuable to all the stakeholders in the VR ecosystem. As a standalone application, this tool can in fact be used to assess the acceptability of the content (both videogames and 360° videos) in the early stages of the production chain or can be used to validate the efficiency of sickness reduction strategies (e.g., reduction of the FOV).

Table 1 .

 1 List of extracted physiological features

	BVP features	EDA features
	HR (Heart Rate)	Mean EDA
	AVNN (Average of all NN intervals)	Amplitude: mean, standard deviation,
	SDNN (Standard deviation of all NN	minimum and maximum
	intervals)	Half of recovery: mean, standard
	rMMSD (Root-mean square differences of	deviation, minimum and maximum
	successive of all NN intervals)	IMF(1 to 3): mean, standard deviation,
	nn50 (The number of interval differences of	minimum and maximum
	successive NN intervals greater than 50 ms)
	pNN50 (Percentage of differences between
	adjacent NN intervals that are greater than
	50 ms)
	VLF (Total spectral power of all NN intervals
	between 0.003 and 0.04 Hz)
	LF (Total spectral power of all NN intervals
	between 0.04 and 0.15 Hz)
	HF (Total spectral power of all NN intervals
	between 0.15 and 0.4 Hz)
	LF/HF (Ratio of low to high frequency power)
	Total Power (Total spectral power)

Table 2 .

 2 Results for Random Forest, Gradient Boosting and Support Vector Machine

	Model	RMSE	MAE	R-Squared
	Random Forest	1.52	1.15	.45
	Gradient Boosting	1.55	1.05	.48
	Support Vector Machine 1.57	1.31	.03

https://www.oculus.com/rift/

https://www.vive.com/eu/

http://www.shimmersensing.com/products/shimmer3-wireless-gsr-sensor

An EDA peak is especially characterized by the amplitude (the height of the peak) and the recovery time (time to return to the level EDA before the peak).

RMSE and MAE evaluate the distance between the predicted data and the ground truth (0 corresponds to a perfect prediction). R-Square corresponds to the explained variance (it ranges from 0 to 1, with 1 indicates that the model explains all the variability of the output variable).