SEMI-SUPERVISED LEARNING WITH STRUCTURED KNOWLEDGE FOR BODY HAIR
DETECTION IN PHOTOACOUSTIC IMAGE
Ryo Kikkawa1 , Hiroyuki Sekiguchi2 , Itaru Tsuge2 , Susumu Saito2 and Ryoma Bise1,2
1

Kyushu University, 2 Kyoto University Hospital

ABSTRACT
Photoacoustic (PA) imaging is a promising new imaging technology for non-invasively visualizing blood vessels inside biological tissues. In addition to blood vessels, body hairs arealso visualized in PA imaging, and the body hair signals degrade the visibility of blood vessels. For learning a body hairclassifier, the amount of real training and test data is limited, because PA imaging is a new modality. To address thisproblem, we propose a novel semi-supervised learning (SSL)
method for extracting body hairs. The method effectively
learns the discriminative model from small labeled training
data and small unlabeled test data by introducing prior knowledge, of the orientation similarity among adjacent body hairs,
into SSL. Experimental results using real PA data demonstratethat the proposed approach is effective for extracting bodyhairs as compared with several baseline methods.
1. INTRODUCTION
Photoacoustic (PA) imaging can noninvasively visualize the
3D structure of blood vessels in vivo, which is useful in earlyclinical diagnosis of cancer and many other diseases. In the
PA imaging process, the imaged tissues absorb laser energy
and convert it into thermoelastic expansion, thereafter emitting ultrasonic waves. The 3D structures of the sound sources
can then be reconstructed by sensing the PA waves [1]. In
addition to blood vessels, body hairs also have a PA characteristic. Even with body hairs shaved before imaging, hairroots under skin surfaces are visualized in the PA image, as
shown in Fig. 1. The body hair signals degrade the visibility of blood vessels. In addition, body hair metrics, such as
the number and thickness of hairs are useful for research on
alopecia and cosmetic surgery. In this paper, we propose abody hair detection method for PA imaging.
To extract body hair regions, conventional image processing approaches, such as the combination of thresholding andmorphological operations, have limitations, because the intensity features of blood vessels and body hairs are similar
and they often touch each other along with PA artifacts, as
shown in Fig. 1. Supervised machine learning methods for
3D segmentation, such as V-net [2], usually require a large
number of labeled samples. It takes a high annotation cost,
however, to obtain such large annotation data from 3D vol-Fig. 1. Example of 3D volume captured by PA imaging. Orientations in a local area tend to be oriented in the same direction. Body hairs touch with vessels, where a red circle
indicates the touching area.
ume images, and the amount of data captured from humanbodies is limited, because PA imaging is a new technology.
To effectively learn a model from a small training dataset,semi-supervised learning (SSL) methods have been proposed[3]. To improve classification performance, SSL exploits the
distribution information of large unlabeled data in the feature
space. Because the amount of unlabeled data captured by new
PA devices is also limited, the current SSL approaches are notsufficient to address our problem.
On the other hand, in medical images, prior knowledge
such as the structure of an entire image is often useful for detection and classification. In our case, body hairs under theskin surface in a local area tend to be oriented in the same
direction, as seen in Fig. 1. In conventional SSL, the feature
space is extracted from a single target object but not from a
global spatial structure, such as combination features of only
positive samples. In our case, the orientation distribution differs among different test images, and the positive samples are
unknown in the test data. Therefore, it is difficult to directly
incorporate the orientation features of only positive samples
in an SSL framework, because the candidate objects also include many negative samples and their labels are unknown. It
is not trivial to introduce such combination features of only
positive samples into an SSL framework.
Contribution : Hence, the main contribution of this work is
to propose a novel SSL method for extracting body hairs in
PA imaging. The method effectively learns the discriminative

model from small labeled training data and small unlabeled
test data by introducing prior knowledge (i.e., the orientation
similarity among adjacent body hairs) into SSL. To achieve
this goal, we first investigated the characteristics of body hairregions in PA images and designed object-level features thatare effective for identifying body hairs. In addition, to introduce the structured knowledge that is combination features of
only positive samples, the method alternatively learns the regressor for the feature space extracted from each sample andthat for the feature space expressing the orientation similarity among adjacent body hairs, such that each regressor estimates the likelihood that a candidate region is a body hair.
Our method requires only small labeled data, and each label for a body hair only requires to annotating the two endpoints of the hair. This annotation process is much less costly
than pixel-level annotation in a 3D volume. To the best of
our knowledge, this is the first attempt to automatically extract body hair regions in PA images. Experimental results
using real PA data demonstrate that the proposed approach iseffective for extracting body hairs as compared with several
baseline methods.
Related works : Many SSL methods have been proposed toeffectively use a large amount of unlabeled data [3][8]. For
example, the self-training approach first trains with a small
amount of labeled data, and then confident unlabeled data isadded to the training data and the classifier is trained iteratively [4]. Graph-based SSL expresses the input distribution
by a graph in which the nodes are labeled and unlabeled samples and edges reflect the similarities of examples. Based on
the assumption that similar data points have similar labels,
graph-based SSL then propagates labels from labeled data
nodes to unlabeled data nodes [5][6][7]. These methods usually use only object-level features extracted from each sample
rather than globally structured features. Co-training methods[9] iteratively use the predictions made by models trained on
different views of the same data to label the unlabeled set andupdate the model with the predicted labels. The co-training
approach is based on the principle of maximizing the consensus among multiple independent hypotheses. Brefeld et.al.
[10] introduced this principle into a semi-supervised support
vector learning algorithm for joint input-output spaces in thefield of natural language processing (NLP). None of thesemethods consider how to handle features extracted from thespatial structure, such as combination features of only positive samples, in the SSL framework.
2. SSL USING STRUCTURE KNOWLEDGE
Figure 2 shows an overview of our SSL method. The proposed method first detects a set of redundant candidate regions. Then, it identifies whether a candidate region is a bodyhair. To effectively use prior knowledge (the orientation similarity among adjacent body hairs), the method alternately
trains a regressor for the object-level features and a regressor for the features expressing the distribution of orientation

Training
Image

Test
Image

Candidate
RegionDetection
Positive
Negative
Unlabeled
Intensity-BasedFeature Space 𝐹

LabeledSamples 𝐴

Unlabeled
Samples 𝑋

Z = 𝑋𝑗 , 𝑞(𝑋𝑗 ; Θ𝑓 , Θ𝑔 )Train
Regressor𝑓(∙ ; Θ𝑓 )𝑋𝑗 , 𝑓(𝑋𝑗 ; Θ𝑓 )𝑤𝑓 𝑓 𝑥; Θ𝑓
Train
𝑤𝑔 𝑔 𝑥; Θ𝑔Regressor𝑔(∙ ; Θ𝑔 )Joint
Regressor𝑞(∙ ; Θ𝑓 , Θ𝑔 )CheckConvergence

WeightUnlabeled
Data by 𝑞

Orientation-BasedFeature Space 𝐺

Output𝑋𝑗 , 𝑞(𝑋𝑗 ; Θ𝑓 , Θ𝑔 )Fig. 2. Illustration of the overview of our SSL method.
similarity among adjacent body hairs. In each iteration, unlabeled data is weighted by the current regressors, and both thelabeled training data and the weighted unlabeled data are used
for updating both regressors. The following sections explain
the details of each step.
2.1. Candidate region detection
The goal in this step is to produce a set of candidate regionsthat may include many false positives but very few false negatives. This would indicate that most true positives are included in the candidate set. As shown in Fig. 1, body hairs
often touch blood vessels in PA imaging, because of radial artifacts in which the intensities of the regions between the bodyhairs and blood vessels are higher than the background intensity. Because the intensity of an artifact is usually slightly
less than that of a sound source (i.e., a vessel or hair), candidate regions are identified by segmenting all regions throughmultiple-level thresholding. We set K thresholds that are
equally spaced, and each threshold is used to segment images
at a particular level of intensity. The detected regions constitute a tree structure with the candidates as nodes. To reduce
the number of candidates, if a node has only one child, thechild node is pruned [11]. In the training data set, if a detected candidate region includes two annotation points placedat the endpoints of a body hair, it is labeled as a positive sample; otherwise, it is labeled as negative. We denote the set ofllabeled candidate regions as A = {Ai , li }N
i=1 , where Nl isthe number of regions, and li is the i-th region’s label. In atest image, all candidate regions are used for unlabeled dataudenoted by X = {Xi }N
i=1 .

2.2. Regressor in object-level feature space
In this section, we explain the object-level features that are
extracted from each candidate region. To design the objectlevel features to enable identification of body hair or otherregions, we first investigated the characteristics of PA imaging. In general, a PA spectrum of a sample can be recordedby using different wavelengths of light. This spectrum canbe used to identify the absorbing components of the sample.
In such multi-spectrum images, we can observe the following characteristics: (1) the intensities of a body hair region
tend to be larger than those of a blood vessel region; (2) the
difference between signals captured using two different wavelengths of light (e.g., 756 and 797nm) tends to be larger in abody hair region than in other regions; and (3) the major axis
length tends to be larger for a body hair region than for otherregions.
Given these characteristics, we use the intensity histogram, the intensity gradient histogram, the intensity histogram of the difference between images captured by two
different wavelengths, and the major axis length of the candidate region as object-level features. Then, we train the
regressor f ( · ; Θf ) in the object-level feature space F byusing logistic regression, in which the regressor estimates
the likelihood that a candidate region is a body hair. The
regressor trained with the labeled data A is used as an initial regressor in the subsequent step. In the inference for atest image, this regressor is iteratively trained and updated
using the labeled data and weighted unlabeled data in the test
image.
2.3. Regressor in structured-prior-based feature space
We here define the feature space G that expresses the orientation similarity among adjacent body hairs, and the regressor g
for this space. To express the orientation similarity in a local
area, we use a sliding window method to separate the entire
volume into M local regions denoted as L = {L1 , ..., LM },where the window size is 1/16 of the entire volume, and thesliding distance is half the window size. This indicates that
G and g consist of M spaces, G = {G1 , ..., GM }, and M
regressors, g = {g1 , ..., gM }, respectively. In the test image,when the center of Lm is the nearest to Xj , each unlabeled
candidate Xj is assigned to the corresponding space Gm .
We next explain each feature space Gi and the corresponding regressor gi . The feature space Gi expresses theorientation distribution of candidate regions in the corresponding local window Li . To compute the orientation of thej-th candidate region, we use the normalized first principal
component (xj , yj , zj ) of the region as the orientation vector,
and then all orientation vectors in Li are mapped to space
Gi , as illustrated in Fig. 3. Because the orientation vector is
normalized to have a norm of 1, the points are distributed on
the surface of a sphere of radius 1. To measure distance on
the sphere’s surface, we map the 3D space to a 2D space byusing stereographic projection as shown in the figure. In this
distribution, because the orientations among adjacent body

Window Size
280 × 280 ×400

𝐺1

𝐿1Sliding Window
by 140 Pixels

𝐺2

𝐿2

1120

⋮
1120

Fig. 3. Examples of local windows and their corresponding
feature spaces, G. Left column: maximum intensity projection (MIP) of entire volume. Middle column: MIPs of local
windows. Right column: normalized orientation distributions
and their stereographic projection.
hairs tend to be similar, we assume that the standard deviationof positive samples (body hairs) becomes smaller than thatof negative samples (other regions), given approximation by
a Gaussian distribution. From this assumption, if the labelof each sample is known, then we can simply estimate theposterior P (hair|Xj ) with a Bayesian classifier:g(Xj ; Θg ) = P (hair|Xj )=p(Xj |hair)P (hair)
P (hair)p(Xj |hair) + P (other)p(Xj |other)(1)

where p(Xj |hair) and p(Xj |other) are estimated fromGaussian distributions of the positive and negative samples,
respectively, and P (hair) and P (other) are the proportionsof the samples. Θg is a set of parameters for modeling the
distributions.
In the test step, however, the label of each sample is unknown. Instead of labels, we thus use the likelihood obtained
from the object-level feature regressor f in the first iteration
step. Each candidate region is weighted using f (Xj ), andthen the Bayesian classifier is iteratively updated with theweighted unlabeled data. Note here that the orientation distribution of body hairs in the training (labeled) data may differ
from that in the test (unlabeled) data. Thus, we only use theweighted unlabeled data in the inference for a test image.
2.4. Alternative semi-supervised-learning in inference
Here, we give an overview of our method. In the pre-training
step, the method trains the object-level regressor f ( · ; Θf ) byusing the labeled data A. In the inference step, the method
first detects candidate regions as unlabeled samples X fromthe test data. Given A and X, it infers the likelihood thatan unlabeled sample is positive (body hair), for all unlabeled
samples X by using the regressor f ( · ; Θf ) pre-trained withA, and the likelihood f (Xj · ; Θf ) is used for the weights toTable 1. Accuracy of compared methods.
Precision Recall F1 score
Supervised
0.908
0.958
0.930
Semi-supervised (0.9)
0.915
0.961
0.935
Semi-supervised (0.6)
0.917
0.966
0.940
Proposed
0.923
0.985
0.952

Fig. 4. Examples of virtual hair removal results, with theimage contrast manually adjusted for visualization purposes.
Left column: original PA images. Right column: results after
hair removal.
train g( · ; Θg ). Then, each unlabeled sample Xj is weighted
with the joint likelihood q(Xj ; Θf , Θg ) = wf f (Xj ; Θf ) +wg g(Xj ; Θg ), where wf and wg are parameters for the user
setting. The set of weighted unlabeled data is denoted asu
Z = {Xi , q(Xj )}N
j=1 . Using both the labeled data A andthe weighted unlabeled data Z, the method re-trains the regressors f ( · ; Θf ) and g( · ; Θg ). This procedure is iterated
done until it converges or reaches maximum number of iterations. Empirically, it converged in less than ten iterations inall our experiments.
3. EXPERIMENTAL RESULTS
We first quantitatively evaluated our method on a real 3D PA
volume dataset captured from multiple patients, in which eachvolume had dimensions of 1120 × 1120 × 400. To generate
the ground truth, we annotated the endpoints of body hairs
for three real 3D PA volumes. We set the parameter K as 20,
which did not affect the result but just changed the computational cost. Then, the candidate regions and annotated points
were matched to generate labeled training data. In all experiments, we simply set wf = wg = 0.5.
Because there is no existing method for automatically extracting body hair regions in PA images, we also developed
two baseline methods to show the effectiveness of the proposed method. These two baselines were also based on thesame candidate detection and classification approach as the
proposed method. Baseline 1 was a supervised method thattrained the logistic regressor in the proposed intensity-based
feature space by using only the training data A; it then classified the unlabeled samples by using the trained regressor.
Baseline 2 was a semi-supervised method that first pre-trained
the regressor by using the training data and then iteratively
self-trained the regressor by adding confident data to the train-

ing data in the inference step. As this method has the confidence threshold as a user-set parameter, we applied several
thresholds (0.6, 0.7, 0.8, 0.9) and then selected the best case(0.6) and worst case (0.9) for comparison.
To evaluate these methods, we performed cross-validationleaving one volume out and computed the average precision,
recall, and F1 score. Table 1 lists the results. In this evaluation, the proposed method achieved the best performance by
all metrics. The results show the effectiveness of the proposed
method for training the classifier from small labeled training
data and small unlabeled test data.
In addition to the above quantitative evaluation, we applied our method to virtual hair removal to improve the visibility of blood vessels in PA imaging. We first classified body
IPh17Lhair regions by the proposed method. Then, if a candidate region was classified as a body hair, we filled the region withthe average intensity of the periphery of the region. This procedure was applied to eight real PA images that were unlabeled. Figure 4 shows examples of the virtual hair removal
results. The left column shows maximum intensity projection(MIP) images of the original 3D PA volume. The right column shows the resulting images after body hair removal. Notethat the high-intensity points are not body hairs but markers for identifying position, as shown by the red circle in theupper-right image. These points were correctly classified as
negatives by our method. In these experiments, the proposed
method worked well for removing body hairs, even though thetest data included diverse conditions, such as different numbers of body hairs and different intensity distributions.
4. CONCLUSION
In this paper, we proposed a novel SSL method for extracting body hair regions in PA imaging. The method can effectively learn a discriminative model from small labeled training data and small unlabeled test data by introducing priorknowledge, consisting of the orientation similarity among adjacent body hairs, into SSL. Experimental results using real
PA data demonstrated that the proposed approach is effective
for extracting body hairs as compared with several baseline
methods.
Acknowledgments.
This work was supported by ImPACT Program of Council for Science, Technology and Innovation (Cabinet Office,
Government of Japan) and JSPS KAKENHI Grant Number
JP18025450.

5. REFERENCES[1] Li, C., and Wang, L.V.; Photoacoustic tomography andsensing in biomedicine, Phys Med. Biol. 54 (19), pp. 5997 (2009).
[2] Milletari, F., Navab, N., and Ahmadi, S.A.; V-Net: Fully
convolutional neural networks for volumetric medical image segmentation, 3DV, pp. 565-571 (2016)[3] Zhu, X.; Semi-supervised learning literature survey,
Computer Sciences TR 1530, University of Wisconsin Madison (2007).
[4] Rosenberg, C., Hebert, M., and Schneiderman, H.; Semisupervised self-training of object detection models, IEEE
WACV, vol. 1, pp. 29-36 (2005).
[5] Blum, A., Lafferty, J., Rwebangira, M. R., and Reddy,
R.; Semi-supervised learning using randomized mincuts,
ICML, pp. 13-20 (2004).
[6] Subramanya, A., Petrov, S., and Pereira, F.; Efficientgraph-based semi-supervised learning of structured tagging models, in EMNLP, pp. 167-176 (2010).
[7] Kingma, D.P., Rezende, D.J., Mohamed, M., andWelling, M.; Semi-supervised learning with deep generative models, NIPS, pp. 1-9 (2014).
[8] Li, Y., and Zemel, R.; High order regularization for semisupervised learning of structured output problems, ICML
vol. 32, pp. 1368-1376 (2014).
[9] Blum, A., and Mitchell, T.; Combining labeled and unlabeled data with co-training, COLT (1998).
[10] Brefeld, U., and Scheffer, T.; Semi-supervised learningfor structured output variables, ICML (2006).
[11] Bise, R., and Sato, Y.; Cell detection from redundantcandidate regions under non-overlapping constraints,
IEEE TMI, vol. 34, pp. 1417-1427 (2015).