Improving the measurement of the Higgs boson-gluon coupling using convolutional neural networks at $e^+e^-$ colliders
Gexing Li, Zhao Li, Yan Wang, Yefan Wang

TL;DR
This paper demonstrates that convolutional neural networks significantly improve the precision of measuring the Higgs boson-gluon coupling at lepton colliders by analyzing jet constituents, reducing uncertainty by approximately 35%.
Contribution
The study introduces CNN-based techniques for Higgs-gluon coupling measurement, outperforming traditional methods in collider data analysis.
Findings
Uncertainty reduced from 1.94% to 1.28% with PYTHIA data.
Uncertainty reduced from 1.82% to 1.22% with HERWIG data.
CNNs utilizing jet constituent energy distributions enhance identification accuracy.
Abstract
In this paper we propose to use convolutional neural networks (CNNs) to improve the precision measurement of the Higgs boson-gluon effective coupling at lepton colliders. The CNN is employed to recognize the Higgs boson and a boson associated production process, with the Higgs boson decaying to a gluon pair and the boson decaying to a lepton pair at the center-of-mass energy 250 GeV and integrated luminosity 5 ab. By using CNNs, the uncertainty of the effective coupling measurement can be decreased from to about using the PYTHIA data and from to about using the HERWIG data in the Monte Carlo simulation. Moreover, the performance of CNNs using different final state constituents shows that the energy distributions of the leading and subleading jets constituents play a major role in the identification and the optimal uncertainty ofâŠ
| PYTHIA | 2.93% | 2.53% | 1.94% |
| HERWIG | 2.67% | 2.47% | 1.82% |
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Improving the measurement of the Higgs boson-gluon coupling using convolutional neural networks at colliders
Gexing Li
Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
School of Physics Sciences, University of Chinese Academy of Sciences, Beijing 100039, China
ââ
Zhao Li
Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
School of Physics Sciences, University of Chinese Academy of Sciences, Beijing 100039, China
ââ
Yan Wang
College of Physics and Electronic Information, Inner Mongolia Normal University, Hohhot 010022, China
Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
ââ
Yefan Wang
Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
School of Physics Sciences, University of Chinese Academy of Sciences, Beijing 100039, China
Abstract
In this paper we propose to use convolutional neural networks (CNNs) to improve the precision measurement of the Higgs boson-gluon effective coupling at lepton colliders. The CNN is employed to recognize the Higgs boson and a boson associated production process, with the Higgs boson decaying to a gluon pair and the boson decaying to a lepton pair at the center-of-mass energy 250 GeV and integrated luminosity 5 ab*-1*. By using CNNs, the uncertainty of the effective coupling measurement can be decreased from to about using the PYTHIA data and from to about using the HERWIG data in the Monte Carlo simulation. Moreover, the performance of CNNs using different final state constituents shows that the energy distributions of the leading and subleading jets constituents play a major role in the identification and the optimal uncertainty of effective coupling using CNNs is reduced by about compared to that using conventional method.
pacs:
07.05.Mh, 14.70.Dj, 14.80.Bn
I Introduction
The Higgs boson occupies a distinct place in the Standard Model (SM) of particle physics. Many lingering physics problems are linked to the Higgs boson, for instance, the stability of the vacuum, electroweak hierarchy problem and dark matter. These problems imply the existence of new physics beyond the SM and require a good understanding of the Higgs properties. The effective coupling of the Higgs boson to a gluon pair is one of the most important parameters. Many theories beyond the SM predict that the Higgs boson-gluon coupling may have deviation from the SM prediction by direct or indirect effects, for example, the stop in supersymmetry or the quark in little Higgs models can contribute to the coupling through the loop effects Blum et al. (2013); Han et al. (2013); Einhorn (1993); Kanemura et al. (2004); He et al. (2013); Moyotl et al. (2016); Baek and Yuan (2017); Hou and Kikuchi (2017); Kanemura et al. (2017); PaĂehr and Weiglein (2017). Therefore, the precision measurement of the Higgs boson-gluon coupling will be a touchstone of the SM and may lead to a breakthrough for new physics.
Although the gluon fusion is the most important process of the Higgs boson production at the CERN Large Hadron Collider, the Higgs boson-gluon coupling is still difficult to be determined accurately due to the overwhelming large QCD radiation Peskin (2012, 2013). The better candidates for the precision measurement of Higgs boson-gluon coupling can be electron positron colliders, which have the clean environment and the high luminosity. The possible future electron positron colliders, which are usually called the Higgs factory at 250 GeV center-of-mass energy, include the Circular Electron-Positron Collider CEPC Study Group (2018a, b); Mo et al. (2016), Future Circular Collider-electron-positron Bicer et al. (2014); Barletta et al. (2014); Benedikt and Zimmermann (2016) and International Linear Collider Behnke et al. (2013); Baer et al. (2013); Adolphsen et al. (2013a, b); Abramowicz et al. (2013). At the Higgs factory, the measurement on most of the Higgs properties can reach percent level accuracy Peskin (2012, 2013); Gao (2018). For the Higgs boson-gluon effective coupling the He et al. (2013); CEPC Study Group (2018b) is always used to parameterize its deviation from the SM prediction, where . With the conventional method (only using the kinematic cuts and tagging) Bai (2017) the uncertainty of the will reach about for the channel of a boson decaying to a lepton pair including the detector effect at the Circular Electron-Positron Collider.
The measurement accuracy of the Higgs boson-gluon coupling can be further improved through an effective identification of jet types. In the last few decades, many different observables motivated by color charge, color connections, electrical charge, or spin have been proposed and achieved good performance Gallicchio and Schwartz (2013); Shelton (2013); Larkoski et al. (2017). For example, the jet energy profile is one of the useful jet substructure observables to distinguish quark and gluon jets by the energy distribution of jet constituents. By using the jet energy profile, the uncertainty of the Higgs boson-gluon coupling can be further reduced to about for the channel of a boson decaying to a lepton pair Li et al. (2018).
However, an observable usually only describes a certain aspect of the jets or some special processes. Although it is better to choose a set of complementary observables to extract more comprehensive characteristics to identify different types of jets or events, the applicable scope of different observables and the degree of association between them will also be difficult problems. Moreover, the deeper correlations between the jet or event constituents may be difficult to be extracted by the artificial observables.
Deep learning has been applied to solve many complicated problems in particle physics. In particular, deep neural networks have been employed to distinguish different types of jets, including Higgs boson tagging Chiappetta et al. (1994), boosted boson tagging Cogan et al. (2015); de Oliveira et al. (2016), boosted top tagging Almeida et al. (2015); Pearkes et al. (2017), single merged jet tagging Baldi et al. (2016), heavy-light quark discrimination Guest et al. (2016) and quark-gluon discrimination Lonnblad et al. (1990, 1991); Peterson et al. (1994); Komiske et al. (2017). They all get an exciting recognition capability and superior to the conventional method. A convolutional neural network (CNN) is one of the most popular and powerful algorithms. Its powerful ability of image recognition makes it easy to extract more comprehensive and deeper features to analyze the jet substructure. It is very suitable for jet tagging and also for testing different shower and hadronization schemes by comparing different Monte Carlo (MC) generators.
In this paper, we propose to use the CNN for the precision measurement of Higgs boson-gluon effective coupling by distinguishing the background processes from the process of a boson decaying to a lepton pair and a Higgs boson decaying to a gluon pair () at lepton colliders. The global information in an event is used for the training of the CNN instead of the jet information. We will use events from different event generators for neural network training and testing to illuminate the difference between the different shower and hadronization schemes.
The content is organized as follows. In the next section, the CNN is briefly reviewed. In the third section, the MC events are generated by PYTHIA and HERWIG. The production of images and CNN architecture are introduced in the fourth section. In the fifth section, we show the results using the CNN. The conclusion is made in the last section.
II Convolutional Neural Networks
A neural network is one of the most popular algorithms in machine learning. Generally, a neural network consists of an input layer, hidden layer, and output layer. A layer is dense if each of its units connects to all of the units in the previous layer. If a neural network consists of a dense layer completely, it will tune a large number of parameters and waste a lot of computing resources. Actually, each neuron only needs to perceive the local image instead of the global image for image recognition, and then the global information can be obtained by integrating the local information at a higher level. This motivates the design of the CNN Lecun et al. (1989). In the last few years, based on the development of computer technology, the CNN has been a mainstay of many major breakthroughs in various fields.
In the image identification, the images in the CNN will pass a convolutional layer, pooling layer, and dense layer. The function of the convolutional layer is extracting features of the image. This can be implemented by the convolution of the filter and the image. A filter is a grid of weights, where is the filter size. The convolution is that each weight in a filter multiplies the corresponding pixel intensity in a patch the same size as an image. Then, we sum the convolutional values, add a bias, and feed it to an activation function. Activation functions introduce the nonlinear properties into neural networks, which enable the neural networks to learn the deeper information. The most used activation function in CNNs is rectified linear units (ReLU), which is defined as . Each convolutional layer usually has many different filters to extract different features of an image. For the multichannel images, there are different colors and convolutional filters in each channel. Each color or channel will be solved by a corresponding filter, like the single color image, and will be accumulated in the final step.
Then, a pooling layer, following the convolutional layer, is used to reduce the number of parameters. The filter of the pooling layer is a grid, where is the pooling size. The max pooling and average pooling are the most common pooling functions. Max pooling takes the largest value while average pooling takes the average of all values in a filter region. A dropout usually is added to avoid the overfitting. It refers to the randomly discarding of some neural network units at certain probability in each training Srivastava et al. (2014). Finally, the dense layers are added to integrate the features in the feature maps extracted by the convolution layers and pooling layers to obtain the high-level meanings of the features and then use them for image recognition.
The error of the model can be quantified by the binary cross entropy loss function Chollet et al. (2015)
[TABLE]
where is the number of training events. The and are the real value and the predicted value by the CNN of the th event. The training process is tuning the parameters in the model to minimize the loss function.
III Pre-Processing
The main process of the Higgs boson production is at the future colliders. We choose the process of the boson decaying to a lepton pair and the Higgs boson decaying to a gluon pair () as the signal process since the boson can be reconstructed very well by the lepton pair. The process of different boson decay modes and are discussed first. Then the two lepton channels are combined as . The backgrounds are divided into two-fermion leptonic (final states are a lepton pair from the or intermediate states), two-fermion hadronic (final states are two quarks), four-fermion leptonic (final states are four leptons from the vector boson pair intermediate states), four-fermion semileptonic (final states are a pair of charged leptons and a pair of quarks from the vector boson pair intermediate states), four-fermion hadronic (final states are four quarks), and the Higgs boson production with the final states, which are different from the signal [mainly the Higgs boson and a boson associated production process with the boson decaying to a lepton pair and the Higgs boson decaying to a quark pair () or boson pair ()] Yan et al. (2016); Xin Mo (2017). Both the signal and background events are simulated at future colliders CEPC Study Group (2018a, b); Mo et al. (2016); Bicer et al. (2014); Barletta et al. (2014); Benedikt and Zimmermann (2016); Behnke et al. (2013); Baer et al. (2013); Adolphsen et al. (2013a, b); Abramowicz et al. (2013) for the center-of-mass energy 250 GeV and integrated luminosity 5 ab*-1*. The parton level MC events are generated by WHIZARD 1.95 Kilian et al. (2011); Moretti et al. (2001) and transferred to hadron level by PYTHIA 6 Sjostrand et al. (2006) and HERWIG 7 Bellm et al. (2016), respectively. For clarity, we call them PYTHIA data and HERWIG data, respectively.
We select a pair of isolated leptons to reconstruct the boson. The rest of the final state constituents are clustered into jets via FASTJET 3.3.0 Cacciari et al. (2012) using the anti- algorithm with a large jet cone of , and the energy of each jet is required to be more than 5 GeV. To suppress the two-fermion leptonic and four-fermion leptonic backgrounds Xin Mo (2017), we add two cuts at first. One is the number of the stable charge particles in the final state , and another is the electromagnetic energy ratio in the final state . Then, the kinematic cuts, i.e., invariant mass, recoil mass, and other constraints of the lepton pair and jet pair, are used to ensure that the lepton pair and jet pair, respectively, come from the boson and the Higgs boson to reject the two-fermion hadronic and four-fermion hadronic backgrounds. More details of the analysis can be found in Ref.Li et al. (2018). The reference also shows that the tagging cannot decrease the uncertainty effectively since its mistag rate for the gluon jet will exclude some gluon jets. Therefore, we only use the tagging in this paper.
The kinematic cuts and tagging can remove a large number of the distinct backgrounds, which will greatly improve the efficiency of the neural network. The remaining backgrounds contain the , , , and four-fermion semileptonic. The jets in the backgrounds and four-fermion semileptonic are mainly heavy quark jets and light quark jets, respectively. But the jets in the backgrounds are jets and light quark jets since quite a few of the light quark jets are merged into the jets with a large jet cone of . It is the complex jet types in the backgrounds that make the signal identification be a challenge.
After all the cuts, the uncertainties of should be evaluated. The evaluation of systematic uncertainties requires a detailed detector study and is unknown yet for the Higgs factory. But the statistical uncertainty of around the SM prediction can be explicitly expressed as
[TABLE]
where and are the numbers of the Higgs boson decaying to gluon pair events and total events, respectively.
In Table 1, the second and the third lines are the uncertainties of with the conventional method using PYTHIA data and HERWIG data, respectively. The difference between the results using PYTHIA data and HERWIG data may come from the different shower and hadronization schemes. The -ordered and the angular-ordered schemes are used for shower effect, and the Lund string and the cluster models are used for the hadronization effect in PYTHIA 6 and HERWIG 7, respectively.
IV Architecture of CNN
For the training of the CNN, we use the combined lepton channel . The entire spherical surface, where the azimuthal angle and the polar angle , is treated as a two dimensional plane image. Each image is designed to have a 66-pixel length in the direction and a 34-pixel length in the direction. The energy of all the final state stable particles is discretized into pixels as our pixel intensity at lepton colliders. The images of the signal process are given the sign one and the other images as the background process are given the sign zero. All the images are divided into the training, validation and test sets in proportion to 8:1:1.
The neural network is implemented by using Keras Chollet et al. (2015) with TensorFlow backend. Our CNN architecture is inspired by the VGGNet Simonyan and Zisserman (2014) architectures and consisted of four iterations of convolutional layers and maxpooling layers shown in Fig.1. Then the feature map is flattened and fed to a dense layer with 128 units. Finally, a dense layer with one unit and a sigmoid activation is added to classify the signal and background processes. Each convolutional layer consists of 64 or 128 filters with filter size and a ReLU activation. The uniform distribution is used to initialize the filters. The stride length of the convolution is 1. The first convolutional layer is set without padding to weaken the influence of the edge information of the image at the beginning while the others are set with padding to keep all the information of the feature map. Each maxpooling layer performs a down-sampling with a stride length of 2. A dropout layer follows each maxpooling layer and the dense layer to avoid overfitting. All the dropout rates of dropout layers are 0.5 except that the first one is 0.25.
The binary cross entropy is used as the loss function. The optimization of training uses the Adam algorithm Kingma and Ba (2014) and the learning rate is 0.0005. The training is set with batch size 128 and 100 epochs and an early stopping patience of 5. Thus, the training will stop early if the value of the validation loss does not go down 5 times 111Example code is provided on https://github.com/zhaoli-IHEP/Higgs-ML..
The receiver operator characteristic (ROC) curve is usually used to quantify the performance of neural networks. A ROC curve is generated by plotting the true positive rate against the false positive rate. The area under the curve (AUC) is defined to compare the overall performance of the neural networks. In this paper, the true positive rate is the signal process () acceptance efficiency and the false positive rate is the mistag efficiency of the background processes.
Then we test the performance of our neural network and compare it to several different neural networks. Fig.2 shows the background rejection rate as a function of the signal acceptance efficiency for the CNN with different architectures. The lines marked as â3-convâ, âAlexâ and âMiniVGGâ represent the performance of the CNN architectures in the Refs.Komiske et al. (2017); Krizhevsky et al. (2012); Guo et al. (2018), respectively. The green dotted line is the result using the neural network, which contains three iterations of a convolutional layer and a maxpooling layer. The blue dash line is the result using the famous AlexNet, which uses a stack of convolutional layers to increase the nonlinearity of the neural network and bigger filter size to increase the receptive field. So, the performance of the AlexNet has a significant improvement compared to that of the 3-conv. The red dash-dotted line is the result using the neural network, which is inspired by the MiniVGGNet architecture but with a bigger filter size in the first two convolutional layers. More iterations of the convolution layer stack further enhance the nonlinearity of the neural network and lead to improved performance. According to the advantages of the VGGNet, our neural network uses a stack of convolutional layers with filter size instead of a single convolutional layer with a big filter size, which can increase the nonlinearity of the neural network and reduce the number of parameters. The black solid line is the result using our neural network architecture, which is better than other three neural network structures for the identification of our signal and background processes.
V Results
In this section, we will present the improvement on the uncertainty archived by using the CNN.
Fig.3 shows the background rejection rate as a function of the signal acceptance efficiency for our CNN. The area under these curves are the AUC values of the different cases. Both training and testing have been applied to the PYTHIA and HERWIG data. For convenience, The symbol âP(H)+P(H)â is used to represent training with the PYTHIA (HERWIG) data and testing with the PYTHIA (HERWIG) data. It can be found that at around the background rejection rate can reach about , meanwhile, the signal acceptance efficiency could still be acceptable. Furthermore, it can be seen that the AUC value of the âH+Hâ is slightly better than that of the âP+Pâ. More specifically, the curves of the P+P and H+H are very similar at the low signal acceptance efficiency region , but the curve of the H+H is higher than that of the P+P at the high signal acceptance efficiency region . In general, the performance of the P+P and H+H are similar, which indicates that the similar performance of the shower and hadronization schemes in PYTHIA and HERWIG.
The âH+Pâ and âP+Hâ are training and testing with different data as a cross-check to illustrate the universality of the CNN model. It makes sense to compare the performance of the CNN models, which are trained with the different data but tested with the same data. The CNN models are universal if their performance are similar. By comparing the âP+P(H)â to the âH+P(H)â in Fig.3, the performance of the CNN model tested with different data is just slightly worse than that tested with same data in all the signal acceptance efficiency region. It means that our CNN models do not have too much overfitting since they are not overly dependent on the certain data.
The different ratios of the remaining signal and backgrounds can be obtained on the ROC curve in Fig.3. The uncertainty of after using the CNN at each point can be expressed as
[TABLE]
Fig.4 presents the uncertainty of after CNN as a function of the signal acceptance efficiency using the PYTHIA and HERWIG data. At the optimal point , can reach about by using the P+P and by using the H+H. Compared to Table 1, it shows that can be further reduced by for the P+P and for the H+H. The results using the H+H is about smaller than that using the P+P. The small difference of the results may come from the different shower and hadronization schemes in PYTHIA and HERWIG. The results of the cross check are slightly worse than that of the training and testing with the same data. Comparing the P+P to the H+P, the uncertainties of using the H+P is slightly worse than that using the P+P. But the difference of the P+P and the H+P is less than , which far exceeds the measurement accuracy of the future electron positron colliders. The H+H and the P+H are in the same situation. The similar results mean that our CNN models do not have too much overfitting and the results are reliable.
In the previous part, one image is constructed with the information of all the final state stable particles in an event. To gain insight into the improvement by the CNN and find the most important features of the signal and background, different images are constructed with different final state constituents. The following analysis only uses the PYTHIA data. Fig.5 shows the uncertainty of after the CNN as a function of the signal acceptance efficiency using the different images. The line marked as âallâ is the result using the images constructed with the information of all the final state stable particles, and the line marked as âmultijetâ is the result using the images constructed with the information of all the jets clustered by anti- algorithm in an event. The âmultijetâ result is slightly better than the âallâ result in the region . However, the difference of the âallâ and âmultijetâ results is less than at the optimal points and can be ignored. This indicates that the information of jets makes a major contribution to the identification of the signal and background processes. The reason is that most of the information except the jets in an event is the lepton pair, which are very similar in the signal and background processes after using the kinematic cuts. The line marked as âdijetâ is the result using the images only constructed with the information of the leading and subleading jets. The âallâ and âdijetâ results are very similar, which shows that the leading and subleading jets nearly contribute all the features for the CNN. The âmultijetâ and âdijetâ results are also very similar since most of the events only have two jets with a large jet cone of . If the images are constructed only with the leading and subleading jets, the center of the two jets can be chosen as the image center. Then the constituents of the two jets are discretized into pixels to obtain the âdijet translationâ images. By this operation, the jets will not be split into two parts at the margins of the image. It can be seen that the âdijet translationâ and the âdijetâ results are also very similar, which indicates that the symmetry property in the direction has been recognized by the CNN.
After showing that the information of the leading and subleading jets makes a major contribution to the identification of the signal and background processes, we further analyze the contribution of each jet. Fig.6 shows the uncertainty of after the CNN as a function of the signal acceptance efficiency using the different single-jet images. Each single-jet image has the size with the jet cone and is designed to have pixels. The jet axis is chosen at the image center so that there is a complete jet on the single-jet image. The lines marked as âleading jetâ and âsubleading jetâ represent the results using the leading jet images and the subleading jet images, respectively. We can see that the leading and subleading jets are equally important for the identification. Then the leading and subleading jet images as two different channels are combined as the âdijet 2-channelâ by analogy with the recognition of color images, with red, green and blue intensities treated as separate input layers. Compared to the âdijetâ, which puts the leading and subleading jets in one image, the âdijet 2-channelâ removes the relative location information of the two jets. It can be seen that the âdijet 2-channelâ result is just slightly worse than the âdijetâ result, so the relative location information of the jets is not important for this discrimination. From the above analysis, we can conclude that the leading and subleading jets make a major contribution to the identification of the signal and background processes.
In the third section, the analysis shows that the jets in the signal process are mainly gluon jets and the jets in the background processes can be mainly divided into quark jets and jets. The three types of jets have different energy distributions of their constituents. Each jet image using energy as pixel intensity records the energy distribution of the jet constituents. This information can be extracted from the jet images by the CNN to identify the signal and background processes. Therefore, the energy distributions of the leading and subleading jets constituents make a major contribution to the identification of the signal and background processes.
Fig.7 shows the best result using the CNN (the line marked as âmultijetâ) and the result using the conventional method (the line marked as âconventionalâ) for the PYTHIA data. Comparing to the result using the conventional method, the CNN has a significant improvement in a wide signal acceptance efficiency region. At the optimal point , the uncertainty of can be decreased from to about by using the CNN and reduced by about compared to that using the conventional method for the PYTHIA data. Moreover, the result using the HERWIG data is similar to that using the PYTHIA data.
VI Conclusions
In this paper, the CNN is used to improve the precision measurement of the Higgs boson-gluon effective coupling at lepton colliders. By using the CNN the uncertainty of can be decreased from to about using the PYTHIA data and from to about using the HERWIG data in the channel of a boson decaying to a lepton pair in the MC simulation for the center-of-mass energy 250 GeV and integrated luminosity 5 ab*-1*. The difference between the expected uncertainties using the PYTHIA and the HERWIG data is less than . Moreover, the performance of the CNN using different final state constituents is proof that the energy distributions of the leading and subleading jets constituents play a major role on the identification and the optimal uncertainty of using the CNN is reduced by about compared to that using the conventional method.
This work was supported by the National Natural Science Foundation of China under Grant No. 11675185. Y.W. is supported by the China Postdoctoral Science Foundation under Grant No. 2016M601134. The authors want to thank Yu Bai, Gang Li, Qiang Li and Manqi Ruan for helpful discussions.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Blum et al. (2013) K. Blum, R. T. DâAgnolo, and J. Fan, JHEP 01 , 057 (2013), eprint 1206.5303.
- 2Han et al. (2013) X.-F. Han, L. Wang, J. M. Yang, and J. Zhu, Phys. Rev. D 87 , 055004 (2013), eprint 1301.0090.
- 3Einhorn (1993) M. B. Einhorn, in Conference on Unified Symmetry in the Small and in the Large Coral Gables, Florida, January 25-27, 1993 (1993), pp. 407â420, eprint hep-ph/9303323.
- 4Kanemura et al. (2004) S. Kanemura, Y. Okada, E. Senaha, and C. P. Yuan, Phys. Rev. D 70 , 115002 (2004).
- 5He et al. (2013) X.-G. He, Y. Tang, and G. Valencia, Phys. Rev. D 88 , 033005 (2013), eprint 1305.5420.
- 6Moyotl et al. (2016) A. Moyotl, S. Chamorro, H. Castilla-Valdez, and M. A. Pérez (2016), eprint 1610.06299.
- 7Baek and Yuan (2017) S. Baek and X.-B. Yuan, Phys. Lett. B 774 , 662 (2017).
- 8Hou and Kikuchi (2017) W.-S. Hou and M. Kikuchi, Phys. Rev. D 96 , 015033 (2017).
