# Evaluating Deep Taylor Decomposition for Reliability Assessment in the Wild

Stephanie Brandl, Daniel Hershcovich, Anders Søgaard

###### Abstract

We argue that we need to evaluate model interpretability methods ’in the
wild’, i.e., in situations where professionals make critical decisions, and
models can potentially assist them. We present an in-the-wild evaluation of
token attribution based on Deep Taylor Decomposition, with professional
journalists performing reliability assessments. We find that using this method
in conjunction with RoBERTa-Large, fine-tuned on the Gossip Corpus, led to
faster and better human decision-making, as well as a more critical attitude
toward news sources among the journalists. We present a comparison of human
and model rationales, as well as a qualitative analysis of the journalists’
experiences with machine-in-the-loop decision making.

## Introduction

Deep neural NLP models increasingly assist humans in making documents easier
to find and analyze. Generally, models are used for either batch processing,
e.g., summarizing or indexing documents, or for online decision support, e.g.,
flagging probable fraud or toxic speech. Deep neural models are often thought
of as black boxes, unable to provide rationales for their decisions, but
recently, many methods have been developed for post-hoc interpretation of deep
neural model predictions (Guidotti et al. 2018; Ancona et al. 2018; Poerner,
Schütze, and Roth 2018; Atanasova et al. 2020). Most such methods provide
rationales in the form of input token attributions. While alternatives exist
(Søgaard 2021), we evaluate attributions below.

Rationales in the form of input token attributions have traditionally been
evaluated automatically or through comparison with gold-standard, human
rationales. Automatic evaluation methods look at the impact of removing (or
keeping only) the tokens that are attributed most relevance or influence,
while comparisons with human rationales typically evaluate predicted
rationales in terms of matching metrics such as token-level $F_{1}$.

Figure 1: Celebrity news stories, the text of which (both longer than shown)
are included in the GossipCop corpus. These examples are labeled unreliable.
We ask professional journalists to relabel the texts, with or without the help
of a RoBERTa-based document classifier and Deep Taylor Decomposition
rationales.

### Related work

Horne et al. (2019) evaluated if machine learning with and without feature
attribution rationales can help lay people assess the reliability of claims.
We follow their data collection set-up in most respects, but not their
experimental protocol. We also rely on a different source of data. Horne et
al. (2019) found a positive effect of model rationales on human decision
making, but evaluated this only in terms of mean ratings for reliable and
unrealiable articles, leaving it open whether human decision-making was more
accurate at the article level. They found no positive effect of showing just
model confidences. The most important difference between their study and ours,
however, is that we evaluate the usefulness of reliability detection with
feature attribution rationales on professional experts (journalists).

More recently, Mohseni et al. (2021) evaluated the effect of model rationales
(using self-explanatory models) across four fake news detection models, also
with lay people. They found that using rationales is not significantly better
than showing only model confidence, and generally found the effectiveness of
explanations to vary across architectures. Our results align better with Horne
et al. (2019) than with Mohseni et al. (2021), but provide a more nuanced
picture and for a population of end users, who are trained to do reliability
assessment, and who have a real need to perform this activity on a day-to-day
basis.

In older work, Feng and Boyd-Graber (2019) evaluated the usefulness of model
rationales in the context of human question answering. They recruited trivia
experts and novices to play QuizBowl with a computer as their teammate, with
models returning both top-$k$ predictions and rationales. Rationales was
provided both in terms of feature attributions and training examples. Their
work is also limited in several respects, relying on an inherently
interpretable linear question answering model and therefore, in effect, not
evaluating interpretability methods, but simply the usefulness of an
interpretable model. Their task is also different than ours in having an open-
ended output space. This means that simply listing top-10 candidate answers
can contribute to human recollection of answers, an effect they exploit, but
which we cannot for a binary classification task. Feng and Boyd-Graber (2019)
find that only novices generally benefit from model rationales. Finally,
QuizBowl is arguably not a task where interpretability is needed – in the
sense that this is a necessary feature for decision support in reliability
assessment to avoid the spread of misinformation. In sum, their work leaves
open a) whether interpretability methods are useful when applied to deep
neural models, b) whether interpretability methods can be useful to experts,
and c) whether interpretability methods can be useful in situations where
interpretability is critical. While Horne et al. (2019) and Mohseni et al.
(2021) addressed a) and b), our work provides partial answers to all these
three questions.

Other work that addresses human-in-the-loop evaluation of interpretability for
deep neural models (a) includes Gonzalez and Søgaard (2020) and González et
al. (2021), but both evaluate interpretability methods with lay people and on
non-critical tasks, ignoring (b) and (c). Attempts to evaluate
interpretability methods for experts performing critical tasks, have, to the
best of our knowledge, been limited to automatic evaluation or evaluation
against gold-standard human rationales.

### Contributions

In sum, our main contribution compared to previous work is that we evaluate
feature attribution rationales for reliability assessment on professional
journalists. While misinformation spread in social media is devastating,
misinformation accelerates when authoritative platforms, including traditional
media, pass on unreliable stories (Soares and Recuero 2021). As a way of
evaluating interpretability methods, professional journalists present a much
higher bar. Journalists are already trained to assess the reliability of news
stories, and predictions and rationales have to be of very high quality to
assist, rather than bias, such professional end users. In addition, we also
evaluate a novel, state-of-the-art interpretability algorithm (Deep Taylor
Decomposition for deep neural networks) on a new source of data, namely a
corpus of celebrity gossip stories, a domain which is arguably a greater
source of misinformation than any other domain and has been estimated to
account for up to half of all online misinformation (Acerbi 2019). Finally, we
consider the impact of the demographics and experience of the professional
journalists on the usefulness of feature attribution.

Figure 2: POS distribution of most highlighted words Figure 3: The interface used in our human reliability assessment with professional journalists receiving model prediction (with confidence) and rationales as feedback. | Time | Error
---|---|---
Text | 45.8 | 0.31
\+ Confidence | 49.2 | 0.25
\+ Confidence + Feat.attr. | 48.6 | 0.22
Table 1: Deep Taylor Decomposition (last row) enables faster (Time) and better
(Error) reliability assessment, even for professional journalists.

Figure 4: Error rate for all 3 conditions with journalists grouped by
experience (top), gender (middle), and by how useful they rated the model
feedback (right). Error bars are standard deviation across bootstrapped
samples.

## Experimental Setup

### Participants

We recruited 25 professional journalists (6 female, 17 male, 2 without
response) to participate in our experiment. All journalists were employed at
the time of the experiment in the same media house, working for a tabloid
newspaper, and had in this capacity daily contact with reliability assessment
of celebrity-related news. The participants were compensated through salary.

### Data

We sampled 80 articles for our experiments from the GossipCop section of the
FakeNewsNet (Shu et al. 2019).111github.com/KaiDMML/FakeNewsNet The articles
were split in 2 batches of 40 articles, half of which our model (see below)
had classified incorrectly. In each condition, each article was seen by a
minimum of three journalists. Gossip as a genre has the advantage that
reliability assessment here should be relatively unbiased by domain experience
and political alignment, evaluating instead the plausibility of a reported
state of affairs. We provide an example in Table 1.

### Model

We train a binary document classifier by fine-tuning RoBERTa-Large (Liu et al.
2019) on the rest of the GossipCop data. RoBERTa-Large – a pre-trained English
language model based on Transformer blocks (Vaswani et al. 2017) – is
augmented with a simple classification architecture to calculate cross-
entropy. Our architecture is completely standard and applies drop-out to the
Transformer’s last layer hidden state. The regularized output is then fed into
a single fully-connected layer and a softmax activation function to scale
logits to probabilities. We use Deep Taylor Decomposition (Chefer, Gur, and
Wolf 2021) to assign feature attribution scores to input words. See code for
details. Deep Taylor Decomposition is non-trivial in Transformer models
because propagation involves attention layers and skip connections, but
Chefer, Gur, and Wolf (2021) present a novel, consistent method for computing
attribution scores across such layers and connections.

### Predictive words

Deep Taylor Decomposition assigns feature attribution scores to input words.
Reviewing aggregate statistics across the examples used in our experiments, we
see the POS distribution for the most highlighted words in Figure 2. We see
our model is very sensitive to punctuation, e.g., misplaced commas, and nouns,
e.g., first names without last names. Several of the journalists, when asked
for signs of unreliability, mentioned poor grammar as the main indicator
outside of metadata.

### Stimuli presentation

The feature attribution scores were displayed to the journalists in one of
three conditions; see Figure 3 to see the user interface. The interfaces for
the other conditions were identical, except without color highlights and/or
model predictions and confidence scores.

### Metrics

We evaluate the time used by professional journalists to assess the
reliability of the news articles in terms of average wall-clock time in
seconds. We evaluate their ability to assess reliability by evaluating their
accuracy (1-0 loss) compared to the GossipCop gold-standard annotations.

## Results

Our main results are listed in Table 1. We see that the time journalists take
to assess reliability is not significantly impacted by the model feedback. We
see a 3.5s increase on average when introducing model predictions and
confidence scores, but on the other hand, the highlighting seems to increase
reading speed a bit, leading to 0.6s faster reading times than in the setup
without feature attributions.

More importantly, human error rate is substantially lower for professional
journalists when they receive model feedback. In the baseline condition of no
model feedback, they are only able to correctly estimate the reliability of
two in three documents, whereas with model confidence scores, their accuracy
is three in four, i.e., the error rate is 0.25. Model rationales in the form
of feature attribution further decreases the error rate to 0.22.

## Analysis

As part of our experiment, we also asked the journalists to volunteer
information about their seniority, gender, and how useful they found the model
feedback. We use this additional survey data to perform a more fine-grained
analysis of the usefulness of the feature attribution scores of Deep Taylor
Decomposition for the task of reliability assessment.

### Experience

The 25 professional journalists in our experiment all work for the same media
house, but some are more senior than others. The usefulness of feature
attribution has previously been said to correlate with level of experience
(Feng and Boyd-Graber 2019). Our data corroborates this claim (see Figure 4,
left chart). In fact, it turns out the positive effect of interpretability is
entirely with journalists with less than three years of experience. Model
feedback, in contrast, seems to hurt the reliability assessments of more
experienced journalists.

### Demographics

We also considered whether reported gender of the participants had an effect
on the usefulness of model feedback, i.e., are men or women more prone to
model suggestions? The positive effect of interpretability methods seems
insensitive to the gender of the journalists. See Figure 4 (middle) for
numbers.22223/25 reported they identified as either male or female.

### Usefulness

Finally, we examined whether the perceived usefulness of the model feedback by
the journalists correlated with their error. Mostly, it did not, but the small
group of journalists who found the model feedback very useful had very low
error rates. Most journalists rated model predictions as useful or very
useful, whereas most journalists rated model predictions and rationales as not
useful or rarely useful. This is remarkable in light of the better performance
of journalists relying on both model predictions and rationales. It seems that
showing model predictions and rationales to professional journalists improves
their reliability assessments, but hurts their perception of the usefulness of
having a model in the loop.333This seems paradoxical, but there are many
similar paradoxes in human decision-making, e.g., Simon, Wang, and Keller
(2011).

### Qualitative feedback

The journalists were also asked to provide qualtiative feedback on their
experience with model-supported reliability assessment. Example feedback
ranged from the model and the words marked worked fine for me to semed random
at times. Several journalists explicitly referred to the model as a source of
bias, acknowledging its influence on their assessments.

## Conclusions

We evaluated feature attributions of Deep Taylor Decomposition applied to a
deep neural document classifier based on fine-tuned RoBERTa-Large
representations, in a model-assisted human decision-making experiment with 25
professional journalists performing in-the-wild reliability assessments of
English celebrity news stories. We found feature attribution to have a
positive effect on the ability of journalists to assist the reliability of
these stories. We also saw, however, that this effect was only with
journalists with less than three years of experience. Interestingly, while the
feature attribution scores had a positive net effect on reliability
assessment, a setup where only model predictions and confidence scores were
provided, was perceived of as more useful by the participants. We make the
data and code involved in our experiments publicly available for replication
and future work.444https://github.com/coastalcph/reliability-wild

## Acknowledgements

This work was partially funded by the Platform Intelligence in News project,
which is supported by Innovation Fund Denmark via the Grand Solutions program.
Thanks to Ana Valeria Gonzalez for help with designing our experiments.

## References

  * Acerbi (2019) Acerbi, A. 2019.  Cognitive attraction and online misinformation.  _Palgrave Communications_ , 5: 15.
  * Ancona et al. (2018) Ancona, M.; Ceolini, E.; Öztireli, C.; and Gross, M. 2018.  Towards better understanding of gradient-based attribution methods for Deep Neural Networks.  In _International Conference on Learning Representations_.
  * Atanasova et al. (2020) Atanasova, P.; Simonsen, J. G.; Lioma, C.; and Augenstein, I. 2020.  A Diagnostic Study of Explainability Techniques for Text Classification.  In _Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing (EMNLP)_ , 3256–3274. Online: Association for Computational Linguistics.
  * Chefer, Gur, and Wolf (2021) Chefer, H.; Gur, S.; and Wolf, L. 2021.  Transformer Interpretability Beyond Attention Visualization.  In _CVPR_.
  * Feng and Boyd-Graber (2019) Feng, S.; and Boyd-Graber, J. 2019.  What can ai do for me? evaluating machine learning interpretations in cooperative play.  In _Proceedings of the 24th International Conference on Intelligent User Interfaces_ , 229–239.
  * González et al. (2021) González, A. V.; Bansal, G.; Fan, A.; Mehdad, Y.; Jia, R.; and Iyer, S. 2021\.  Do Explanations Help Users Detect Errors in Open-Domain QA? An Evaluation of Spoken vs. Visual Explanations.  In _Findings of the Association for Computational Linguistics: ACL-IJCNLP 2021_ , 1103–1116. Online: Association for Computational Linguistics.
  * Gonzalez and Søgaard (2020) Gonzalez, A. V.; and Søgaard, A. 2020.  The Reverse Turing Test for Evaluating Interpretability Methods on Unknown Tasks.  In _NeurIPS 2020 Workshop on Human And Model in the Loop Evaluation and Training Strategies_.
  * Guidotti et al. (2018) Guidotti, R.; Monreale, A.; Ruggieri, S.; Turini, F.; Giannotti, F.; and Pedreschi, D. 2018.  A Survey of Methods for Explaining Black Box Models.  _ACM Comput. Surv._ , 51(5).
  * Horne et al. (2019) Horne, B. D.; Nevo, D.; O’Donovan, J.; Cho, J.-H.; and Adalı, S. 2019.  Rating reliability and bias in news articles: Does AI assistance help everyone?  In _Proceedings of the International AAAI Conference on Web and Social Media_ , volume 13, 247–256.
  * Liu et al. (2019) Liu, Y.; Ott, M.; Goyal, N.; Du, J.; Joshi, M.; Chen, D.; Levy, O.; Lewis, M.; Zettlemoyer, L.; and Stoyanov, V. 2019.  RoBERTa: A Robustly Optimized BERT Pretraining Approach.  Cite arxiv:1907.11692.
  * Mohseni et al. (2021) Mohseni, S.; Yang, F.; Pentyala, S.; Du, M.; Liu, Y.; Lupfer, N.; Hu, X.; Ji, S.; and Ragan, E. 2021.  Machine Learning Explanations to Prevent Overtrust in Fake News Detection.  In _Proceedings of the International AAAI Conference on Web and Social Media_ , volume 15, 421–431.
  * Poerner, Schütze, and Roth (2018) Poerner, N.; Schütze, H.; and Roth, B. 2018.  Evaluating neural network explanation methods using hybrid documents and morphosyntactic agreement.  In _Proceedings of the 56th Annual Meeting of the Association for Computational Linguistics (Volume 1: Long Papers)_ , 340–350. Melbourne, Australia: Association for Computational Linguistics.
  * Shu et al. (2019) Shu, K.; Mahudeswaran, D.; Wang, S.; Lee, D.; and Liu, H. 2019.  FakeNewsNet: A Data Repository with News Content, Social Context and Spatialtemporal Information for Studying Fake News on Social Media.  arXiv:1809.01286.
  * Simon, Wang, and Keller (2011) Simon, J.; Wang, Y.; and Keller, L. R. 2011.  _Paradoxes and Violations of Normative Decision Theory_.  American Cancer Society.  ISBN 9780470400531.
  * Soares and Recuero (2021) Soares, F.; and Recuero, R. 2021.  How the mainstream media help to spread disinformation about COVID-19.  _M/C Journal_ , 24(1).
  * Søgaard (2021) Søgaard, A. 2021.  _Explainable Natural Language Processing_.  Synthesis Lectures on Human Language Technologies. Morgan & Claypool Publishers.  ISBN 9781636392141.
  * Vaswani et al. (2017) Vaswani, A.; Shazeer, N.; Parmar, N.; Uszkoreit, J.; Jones, L.; Gomez, A. N.; Kaiser, L. u.; and Polosukhin, I. 2017.  Attention is All you Need.  In Guyon, I.; Luxburg, U. V.; Bengio, S.; Wallach, H.; Fergus, R.; Vishwanathan, S.; and Garnett, R., eds., _Advances in Neural Information Processing Systems_ , volume 30. Curran Associates, Inc.


# Quasi-periodic oscillations of the X-ray burst from the magnetar SGR
J1935+2154 and associated with the fast radio burst FRB 200428

Xiaobo Li Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Mingyu Ge Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Lin Lin  Department of Astronomy, Beijing Normal University, Beijing 100088,
China  Shuang-Nan Zhang* Key Laboratory for Particle Astrophysics, Institute
of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing
100049, China  University of Chinese Academy of Sciences, Chinese Academy of
Sciences, Beijing 100049, China  Liming Song Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Xuelei Cao  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Bing
Zhang  Department of Physics and Astronomy, University of Nevad, Las Vegas, NV
89154, USA  Fangjun Lu Key Laboratory for Particle Astrophysics, Institute of
High Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing
100049, China  Key Laboratory of Stellar and Interstellar Physics and School
of Physics and Optoelectronics, Xiangtan University, Xiangtan 411105, Hunan,
China  Yupeng Xu Key Laboratory for Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049,
China  University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Shaolin Xiong  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Youli Tuo Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Ying Tan
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Weichun
Jiang  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Jinlu Qu  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Shu Zhang  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Lingjun Wang  Key Laboratory for Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049,
China  Jieshuang Wang  Max-Planck-Institut für Kernphysik, Saupfercheckweg 1,
D-69117 Heidelberg, Germany  Binbin Zhang  School of Astronomy and Space
Science, Nanjing University, Nanjing 210023, China  Peng Zhang  Key Laboratory
for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy
of Sciences, 19B Yuquan Road, Beijing 100049, China  Chengkui Li  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Congzhan
Liu Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Tipei Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Department of Astronomy, Tsinghua University, Beijing 100084, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Qingcui Bu  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Ce Cai  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Yong Chen  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Yupeng Chen  Key Laboratory for
Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Zhi Chang  Key Laboratory
for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy
of Sciences, 19B Yuquan Road, Beijing 100049, China  Li Chen  Department of
Astronomy, Beijing Normal University, Beijing 100088, China  Tianxian Chen
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Yibao
Chen  Department of Physics, Tsinghua University, Beijing 100084, China
Weiwei Cui  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Yuanyuan Du  Key Laboratory for Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049,
China  Guanhua Gao  Key Laboratory for Particle Astrophysics, Institute of
High Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing
100049, China  University of Chinese Academy of Sciences, Chinese Academy of
Sciences, Beijing 100049, China  He Gao  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Yudong Gu  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Ju Guan
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Chengcheng Guo  Key Laboratory for Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049,
China  University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Dawei Han  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Yue Huang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Jia Huo  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Shumei Jia  Key Laboratory for
Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Jing Jin  Key Laboratory for
Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Lingda Kong  Key Laboratory
for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy
of Sciences, 19B Yuquan Road, Beijing 100049, China  University of Chinese
Academy of Sciences, Chinese Academy of Sciences, Beijing 100049, China  Bing
Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Gang Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Wei Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Xian Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Xufang Li  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Zhengwei Li  Key Laboratory for Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049,
China  Xiaohua Liang  Key Laboratory for Particle Astrophysics, Institute of
High Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing
100049, China  Jinyuan Liao  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Hexin Liu  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Hongwei Liu  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Xiaojing
Liu  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Xuefeng Lu  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Qi Luo  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Tao Luo  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Binyuan Ma  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  RuiCan Ma  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Xiang Ma  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Bin Meng  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Yi Nang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Jianyin Nie  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Ge Ou
Computing Division, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Xiaoqin Ren  Key Laboratory
for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy
of Sciences, 19B Yuquan Road, Beijing 100049, China  University of Chinese
Academy of Sciences, Chinese Academy of Sciences, Beijing 100049, China  Na
Sai  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Xinying Song  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Liang Sun  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Lian Tao  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Chen Wang  Key Laboratory of Space
Astronomy and Technology, National Astronomical Observatories, Chinese Academy
of Sciences, Beijing 100012, China  University of Chinese Academy of Sciences,
Chinese Academy of Sciences, Beijing 100049, China  Pengju Wang  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Wenshuai Wang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Yusa Wang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Xiangyang Wen  Key Laboratory for
Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Bobing Wu  Key Laboratory
for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy
of Sciences, 19B Yuquan Road, Beijing 100049, China  Baiyang Wu  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Mei Wu  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Shuo Xiao  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Sheng Yang  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Yanji
Yang  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
Qibin Yi  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
School of Physics and Optoelectronics, Xiangtan University, Yuhu District,
Xiangtan, Hunan, 411105, China  Qianqing Yin  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Yuan You  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Wei Yu  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Fan Zhang  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Hongmei Zhang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  Juan Zhang  Key Laboratory for
Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of
Sciences, 19B Yuquan Road, Beijing 100049, China  Wanchang Zhang  Key
Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Wei Zhang
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Yifei Zhang  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Yuanhang Zhang  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China  Haisheng Zhao
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics,
Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China  Xiaofan
Zhao  Key Laboratory for Particle Astrophysics, Institute of High Energy
Physics, Chinese Academy of Sciences, 19B Yuquan Road, Beijing 100049, China
University of Chinese Academy of Sciences, Chinese Academy of Sciences,
Beijing 100049, China  Shijie Zheng  Key Laboratory for Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan
Road, Beijing 100049, China  Dengke Zhou  Key Laboratory for Particle
Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences,
19B Yuquan Road, Beijing 100049, China  University of Chinese Academy of
Sciences, Chinese Academy of Sciences, Beijing 100049, China

###### Abstract

The origin(s) and mechanism(s) of fast radio bursts (FRBs), which are short
radio pulses from cosmological distances, have remained a major puzzle since
their discovery. We report a strong Quasi-Periodic Oscillation (QPO) of
$\sim$40 Hz in the X-ray burst from the magnetar SGR J1935+2154  and
associated with FRB 200428, significantly detected with the Hard X-ray
Modulation Telescope (Insight-HXMT) and also hinted by the Konus-Wind data.
QPOs from magnetar bursts have only been rarely detected; our 3.4 $\sigma$
(p-value is 2.9e-4) detection of the QPO reported here reveals the strongest
QPO signal observed from magnetars (except in some very rare giant flares),
making this X-ray burst unique among magnetar bursts. The two X-ray spikes
coinciding with the two FRB pulses are also among the peaks of the QPO. Our
results suggest that at least some FRBs are related to strong oscillation
processes of neutron stars. We also show that we may overestimate the
significance of the QPO signal and underestimate the errors of QPO parameters
if QPO exists only in a fraction of the time series of a X-ray burst which we
use to calculate the Leahy-normalized periodogram.

methods:statistical – QPO – FRB 200428 – pulsars:individual:SGR J1935+2154  –
X-rays: bursts

##  1 Introduction

A major class of the sources of fast radio bursts (FRBs) (Lorimer et al.,
2007) have long been proposed and recently been proved to be magnetars, which
are neutron stars with surface magnetic field higher than $10^{14}$ Gauss
(Kouveliotou et al., 1998). SGR J1935+2154  is a Galactic magnetar with the
most frequent bursting activities (Israel et al., 2016; Lin et al., 2020b). It
went into burst active episodes in 2014, 2015, 2016, 2019, 2020, 2021, and
2022 (Younes et al., 2017; Lin et al., 2020b, a; Wang et al., 2021; Xiao et
al., 2022). In the active episodes of 2020, a giant radio burst containing two
pulses, which coincides in time with a bright X-ray burst, has been reported
from it (CHIME/FRB Collaboration et al., 2020; Bochenek et al., 2020;
Mereghetti et al., 2020; Li et al., 2021; Ridnaia et al., 2021; Tavani et al.,
2021). The radio burst, denoted as FRB 200428, is the first FRB detected in
other wavelengths, and thus SGR J1935+2154  becomes the first counterpart of a
FRB. This immediately establishes that at least some cosmological FRBs are
produced during magnetar bursts. However, the exact mechanism behind this type
of mysterious phenomena is unclear. It is well established that the majority
of X-ray bursts from SGR J1935+2154  do not come with radio emission down to
very low fluence (Lin et al., 2020c; Kirsten et al., 2021) and weak radio
bursts do not come with bright X-ray bursts (Zhang et al., 2020a; Kirsten et
al., 2021). Therefore this peculiar X-ray burst and its association with FRB
200428 remains the only known FRB-magnetar connection. This X-ray burst is not
a typical magnetar short burst, due to its non-thermal X-ray spectrum (Zhang
et al., 2020b; Li et al., 2021; Younes et al., 2021). Therefore, normal
mechanisms producing magnetar X-ray bursts will not generate (observable)
FRBs, or at least not frequently.

Short X-ray bursts and flares from a magnetar can be generated from starquakes
(Thompson & Duncan, 1995). The damping of crustal oscillations due to core
–crust coupling of the neutron star would leave imprints in the form of Quasi-
periodic Oscillations (QPOs) in the temporal profiles of bursts (Huppenkothen
et al., 2014b; Miller et al., 2019). A series of short-lived QPOs with center
frequencies of $\sim$ 10-600 Hz were found during the giant flares of SGR
1806$-$20 and SGR 1900+14 (Israel et al., 2005; Strohmayer & Watts, 2005,
2006; Huppenkothen et al., 2014b; Miller et al., 2019). Recently, Castro-
Tirado et al. (2021) reported two broad high-frequency QPOs in the main peak
of an extra-galactic giant flare lasting for $3.5\text{\,}\mathrm{ms}$. During
the short bursts from SGR J1550$-$5418, which are in the same category with
those from SGR J1935+2154, a weak QPO signal at $\sim$ 260 Hz was possibly
found in a single burst (Huppenkothen et al., 2014a). Using the same data but
stacking some short bursts together, more significant QPOs centered at
$\sim$93 Hz and $\sim$ 127 Hz, which are close to those found in giant flares,
were reported (Huppenkothen et al., 2014a). Also, a $\sim$57 Hz QPO was found
in stacked short bursts from SGR 1806$-$20 (Huppenkothen et al., 2014c). All
these QPOs can be explained by the oscillations in the crustal movement of
magnetars (Huppenkothen et al., 2014b; Miller et al., 2019). These signals,
though statistically not significant and also not common for most magnetar
X-ray bursts, are among the few links between observations and physics of
neutron stars.

Searching for QPOs in transient light curves requires special cares. Standard
methods involving Fourier analysis are defined for stationary processes.
Magnetar bursts require special cares when performing Fourier analysis on
their light curves because they have a beginning and end in duration and they
are non-stationary processes. For stationary processes, the observed Leahy-
normalized power should follow a $\chi_{2}^{2}$ distribution scaled by the
mean power in each bin, but for the non-stationary signal the distribution at
low frequency (<30 Hz) does not satisfy the $\chi_{2}^{2}$ distribution
(Huppenkothen et al., 2013, 2014a). In this paper, the periodogram is the
squared modulus of the Fourier transform of the light curve. For transient
events where the shape of the burst is known, Monte Carlo simulations of light
curves can be obtained by sampling the fitted profile through a Poisson
distribution (Fox et al., 2001; Guidorzi, 2011). Each light curve can be used
to generate a periodogram, thus one can get the distribution of the
periodogram at different Fourier frequencies. The periodogram distribution at
each frequency is a non-central $\chi_{2}^{2}(\lambda)$ distribution, whose
non-central value $\lambda$ at each frequency is the power spectrum of the
burst profile. This has disadvantages when the burst profiles of magnetars
exhibit a variety of shapes, and thus the profiles could not be determined by
the prior knowledge of the physical origin. This raises the possibility of
false QPO detection or overestimating the significance of the potential QPOs.

There are other methods to search for QPOs in transient time series. A
conservative approach is to model the periodograms by a Bayesian method
(Huppenkothen et al., 2013; Vaughan, 2010). We therefore search for the QPO
signals using the Bayesian method during the short X-ray burst associated with
FRB 200428, observed by the Hard X-ray Modulation Telescope (Insight-HXMT) (Li
et al., 2021) and the light curve of Konus-Wind (Ridnaia et al., 2021). A
significant QPO signal with center frequency at about 40 Hz is detected with
HXMT/ME using the Bayesian method and verified by the Konus-Wind data.

We structure the paper as follows. In section 2, we present the details of the
instrument and the data reduction process for this burst. The dead time
correction of HXMT/ME is introduced as well. In section 3, the result of the
QPO search and its significance is presented using the data of HXMT/ME. The
differences of the significance at different duration of the burst are also
presented and discussed. The cross-correlation analysis between Konus-wind and
HXMT/ME is also presented to verify the QPO. The physics interpretation of the
QPO are discussed in section 4. We give a discussion and summary in section 5.

##  2 The instrument and data reduction

###  2.1 Insight-HXMT /ME

As the first X-ray astronomy satellite of China, Insight-HXMT was launched on
June 15, 2017 and carries three main instruments onboard (Zhang et al.,
2020c): the High Energy X-ray telescope (HE) (Liu et al., 2020), the Medium
Energy X-ray telescope (ME) (Cao et al., 2020), and the Low Energy X-ray
telescope (LE) (Chen et al., 2020). The high flux of the X-ray burst
associated with FRB 200428 caused the data saturation and loss in HE and LE
detectors. The raw data in some time intervals were discarded on-board and
their light curves have several gaps as described in Li et al. (2021).
Therefore, only the data of the ME instrument is used to search for the QPO
and ME also provides the high time resolution (6 us) required for timing
analysis (Tuo et al., 2021).

The data reduction of the X-ray burst associated with FRB 200428 is performed
with the Insight-HXMT Data Analysis Software package (HXMTDAS) version 2.04.
The main steps of ME data reduction are: (1) Use the commands mepical in
HXMTDAS to calibrate the photon events from the 1L data product according to
the Calibration Database (CALDB) of Insight-HXMT. (2) Select the good time
intervals (GTIs) from $T_{0}$ to $T_{0}$+1 s, where $T_{0}$ is 2020-04-28
14:34:24 UTC (Satellite time). (3) Identify the grade of each event according
to their arrival time and calculate the dead time of each FPGA at the
specified interval. (4) Extract the good events based on the GTIs using the
commands mescreen. (5) Produce the light curve of ME with the dead time
correction using the commands, melcgen. Although ME detectors do not have the
problems of the data saturation like HE and LE, the dead time effect is still
significant for bright X-ray burst detection.

We elaborate the dead time effect of ME instrument here. ME consists of three
detector boxes. Each box has three independent Field Programmable Gate Array
(FPGA) and each FPGA manages six Application Specified Integrated Circuit
(ASIC). Each ASIC is responsible for the readout of 32 pixels (Li et al.,
2020). The dead time of ME can be calculated by each FPGA independently. When
one pixel in the FPGA is triggered, a window with $62.5\text{\,}\mathrm{ns}$
is opened and the pixels triggered in other ASICs of the same FPGA will be
recorded. If there are $N_{\rm{hit}}$ of ASICs in the same FPGA that are
triggered within the time window, the electronics will take $163+93\cdot
N_{\rm{hit}}$ microseconds to process these triggered pixels. The typical
timescale of dead time is $256\text{\,}\mathrm{\SIUnitSymbolMicro s}$ because
only one pixel is triggered at most situations (Cao et al., 2020). The dead
time process of ME is similar to that of non-paralyzed detectors, and the
hardware algorithm is implemented in the HXMTDAS to calculate the dead time of
ME. The dead time causes the loss of photons and the distribution of photon
arrivals to deviate from Poisson distribution, and alters the distribution of
power in the periodogram as well. Note that the dead time strongly depends on
the count rate. As a consequence, dead time corrections are especially
important for bright bursts.

A special type of events in ME data could be an indicator to determine whether
the dead time correction is performed well. There are two kinds of thresholds
for ME. One is the triggered threshold, the other is the readout threshold.
The special events are those the pulse height is higher than the triggered
threshold and lower than the readout threshold. These events are mainly the
noise events which are independent of the source brightness, so the light
curve of these events is expected to be constant after dead time correction as
shown in Figure 1.

Figure 1:  Light curve of the special events obtained by Insight-HXMT/ME near
the FRB 200428 burst, where time zero is 2020-04-28 14:34:24 UTC (Satellite
time). The blue line is the light curve after dead time correction and the red
line is the light curve without dead time correction. After dead time
correction, the counts of these special events keep stable.

###  2.2 Konus-Wind

Konus-Wind (KW) consists of two identical NaI(Tl) scintillation detectors,
each with $2\pi$ sr field of view. One detector (S1) points towards the south
ecliptic pole, thereby observing the south ecliptic hemisphere, while the
other (S2) observes the north ecliptic hemisphere. Here the light curves from
KW are obtained from Ridnaia et al. (2021).

###  2.3 Methods for the QPO analysis

####  2.3.1 Monte Carlo simulations of pulse profile

Monte Carlo simulation of light curves is a standard tool in timing analysis.
One can fit an empirical function to the burst profile and then generate a
large number of realizations of that burst profile by adding the Poisson
photon counting noise. The periodograms computed from these simulated light
curves form a distribution against which to compare the periodogram of the
real data. The X-ray burst profile of FRB 200428 are fitted with three
Gaussian functions, which can be found in Li et al. (2021),

$R_{3}=\sum_{i=1}^{3}{N_{{\rm p},i}G(t,t_{{\rm p},i},\sigma_{{\rm p},i})}+l,$
(1)

where $G(t,t_{{\rm p},i},\sigma_{{\rm p},i})=\frac{1}{\sqrt{2\pi}\sigma_{{\rm
p},i}}\exp(-\frac{(t-t_{{\rm p},i})^{2}}{2\sigma_{{\rm p},i}^{2}})$, $N_{{\rm
p},i}$, $t_{{\rm p},i}$ and $\sigma_{{\rm p},i}$ are the normalisation,
arrival time and Gaussian width of the $i$th peak. $l$ is the background level
of the light curve.

Actually, the true pulse profile of the X-ray burst associated with FRB 200428
can not be known precisely due to the degeneracy between the overall burst
profile and the red noise component. If the pulse profile can be presented by
eq. (1), Monte Carlo simulations of pulse profile caused by the Poisson noise
of photon counting can be used to estimate the significance of the QPO
candidate.

####  2.3.2 Bayesian method

In the following, we briefly introduce the Fourier-based periodograms analysed
by Bayesian method commonly used to derive the significance of QPOs (Vaughan,
2010; Huppenkothen et al., 2013). A periodogram is a sample of the power
spectral density (PSD) of the signal based on a given time series
$x_{k}(t)\,(k=0,1,...N-1),N=2^{m}$ (m is an integer). We calculate the
periodogram as the absolute square of the discrete Fourier transform of the
signal, and it follows $\chi_{2}^{2}$ distribution around the PSD
$S_{j}=S(f_{j}),j=0,1,N/2$ at Fourier frequency $f_{j}$ (Van der Klis, 1989;
Timmer & Koenig, 1995; Vaughan, 2010). Therefore, the periodogram $I_{j}$
follows an exponential distribution about $S_{j}$ and can be described below,

$p(I_{j}\mid S_{j})=\frac{1}{S_{j}}e^{-I_{j}/S_{j}}.$ (2)

Many astrophysical transients show excess power at low frequencies, and it is
often assumed that this can be modeled as a red noise process. One basic noise
model $S_{p}(f)$ that is commonly used is a combination of a red noise power
law and a white noise amplitude,

$S_{p}(f)=Af^{-\alpha}+C,$ (3)

where $A$ is the amplitude and $\alpha$ is the index. $C$ is the constant
representing the white noise component. Another noise model $S_{b}(f)$ is
composed of broken power law and a constant white noise,

$S_{b}(f)=Bf^{-\alpha_{1}}+C,(f\leq f_{break})$ (4)

$S_{b}(f)=Bf_{\rm{break}}^{\alpha_{2}-\alpha_{1}}(f)^{-\alpha_{2}}+C,(f>f_{break})$
(5)

where $\alpha_{1}$ and $\alpha_{2}$ are indices of the two power law
components, and $f_{\rm{break}}$ is the break point for the frequency, and $B$
is the normalization value. When $\alpha_{1}$ equals $\alpha_{2}$, the model
of $S_{b}(f)$ will be the same as $S_{p}(f)$ and they satisfy the requirement
of nested models. Following Bayes’ rule, the posterior probability of a set of
model parameters $\theta$, like $A$, $\alpha$ and $C$ in noise model
$S_{p}(f)$, can be defined given the uniform distribution of the prior.

We use the Bayesian method to maximum likelihood estimation (MLE) to the
observed Leahy-normalized periodogram and obtain the maximum a posterior (MAP)
estimates of the model parameters. As introduced in Vaughan (2010), maximizing
the joint likelihood function
$p(\bm{I}|\bm{\theta},H)=\prod\limits_{j}^{N/2}p(I_{j}|S_{j})$ is equivalent
to minimizing the deviance function:

$D(\bm{I},\bm{\theta},H)=-2\log
p(\bm{I}\mid\bm{\theta},H)=2\sum\limits_{j=1}^{N/2}(\frac{I_{j}}{S_{j}}+\log
S_{j}).$ (6)

In other words, finding the maximum likelihood estimates (MLEs) of the
parameters of a model $S(\bm{\theta})$ can be expressed as

$\bm{\hat{\theta}_{\rm{MLE}}}=\arg\min_{\theta}D(\bm{I^{\rm{obs}}},\bm{\theta},H)=\arg\max_{\theta}p(\bm{I}=\bm{I^{\rm{obs}}}\mid\bm{\theta},H),$
(7)

where $\bm{I^{\rm{obs}}}$ is the observed periodogram, $\arg\min$ is to find
the set of parameter $\bm{\theta}$ to obtain the minimum value of function
$D$, and $\arg\max$ of distribution $p$ is to find the maximum likelihood.

To verify the possible QPO in the periodogram, a model selection task is taken
to assess whether the periodogram can be represented by a simple power law or
it requires a more complex model like broken power law. The likelihood ratio
test (LRT) is employed to distinguish two nested model, null hypothesis
$H_{0}$ and the alternative hypothesis $H_{1}$. $H_{0}$ is the power law and
$H_{1}$ is the broken power law. The LRT statistic is twice the minimum log
likelihood of the two models:

$T_{\rm{LRT}}=-2\log\frac{p\left(\bm{I}\mid\bm{\hat{\theta}^{0}_{\rm{MLE}}},H_{0}\right)}{p\left(\bm{I}\mid\bm{\hat{\theta}^{1}_{\rm{MLE}}},H_{1}\right)}=D_{\rm{min}}(H_{0})-D_{\rm{min}}(H_{1}).$
(8)

Once the noise model is selected using the LRT statistic, the QPO search can
also be considered as a model selection problem. The new $H_{1}$ model is the
sum of the selected noise model and a Lorentz model. When the amplitude of the
Lorentz is zero, the $H_{1}$ model becomes the noise model. The asymptotic
theory suggests that the LRT statistic should be distributed as a chi-square
variable under certain regularity conditions,
$T_{\rm{LRT}}\sim\chi_{\nu}^{2}$, where the number of degrees of freedom $\nu$
is the difference between the number of free parameters in $H_{1}$ and
$H_{0}$. When the regularity conditions are not satisfied, one can not use the
distribution to be that of the asymptotic theory (Protassov et al., 2002).
Nevertheless, LRT is still a powerful statistic for comparing models and can
be calibrated by posterior predictive simulations, as shown by Protassov et
al. (2002).

In order to investigate narrow features, another statistic $T_{\rm{R}}$ to use
is the maximum ratio of observed to model power (Vaughan, 2010; Huppenkothen
et al., 2013),

$T_{\rm{R}}={\rm{max}}(R_{j}),$ (9)

where

$R_{j}=2I_{j}/R_{j},$ (10)

$I_{j}$ and $S_{j}$ are the observed and model powers. This is the candidate
for single bin periodicity when searching the periodogram for the highest
data/model outlier. After simulating a large number of periodograms from the
posterior parameter sets and fitting each simulated periodogram with the
preferred noise model, we can also compute the data/model and find the maximum
data/model outlier. Then the distributions of the highest outlier can be
calibrated under the selected noise model. The posterior p-value for the
observed outlier can be obtained when calculating the tail areas from the
simulated outlier distributions.

##  3 Results

After we have obtained the dead time corrected light curve of ME in 18–50 keV
as displayed in Figure 2, the periodograms of this burst is calculated using
the Fast Fourier Transform (FFT) code implemented in open source Python
software stingray (Huppenkothen et al., 2019). We find a potential QPO
candidate with moderate significance at frequency about 40 Hz. We use two
different strategies to evaluate the significance of the potential QPO signal:
(1) Monte Carlo simulation; (2) Bayesian approach. Then the significance at
different duration is also investigated. The cross-correlation between ME and
KW light curves is performed, which also verifies the QPO signal at about 40
Hz.

###  3.1 Monte Carlo simulations of pulse profile

Although Monte Carlo simulation of the pulse profile is an important tool in
timing analysis, it must know the physically motivated burst profile function
in advance. We assume the burst profile of the X-ray burst associated with FRB
200428 can be fitted with three Gaussian functions and the fit parameters can
be found in Table 1. It should be noted that we just take the burst as an
example to show an assumption for the burst profile without taking the red
noise into account will produce an uncorrected significance estimation of the
QPO candidate.

The procedures of the Monte Carlo simulations of the burst profile is
described below. As introduced in section 2.3.1, the burst profile can be
fitted with eq. (1) and shown in Figure 2 with green solid line. The Leahy-
normalized periodograms of light curve can be derived using the Fast Fourier
Transform (FFT) code in Python software stingray as shown in Figure 3. A
potential QPO candidate with moderate significance at frequency about 40 Hz
can be found. In order to establish the p-value of the QPO, Monte-Carlo
simulations of the burst profile are used to validate whether the QPO signal
is generated from the fluctuations of the light curve or not. We then generate
10000 simulated light curves sampled from the fitted three Gaussian burst
profile by adding the Poisson photon counting noise. The periodogram is
computed for each simulated light curve to form a distribution at each
frequency bin. The maximum, mean and minimum values of power distribution at
each frequency bin are calculated and plotted in Figure 3. Comparing the
observed power in each frequency bin with the distribution of simulated powers
allows us to make a statement that the observed power at 40 Hz is difficult to
be produced by the Poisson noise process alone according to the distribution
of simulated powers. From the simulated distribution at each frequency, the
distribution of Leahy power above 10 Hz satisfies a $\chi^{2}_{2}$
distribution. Consequently, the p-value for the QPO at 40 Hz under the
assumption of the three-Gaussian burst profile is 1.7e-7, indicating that the
Poisson noise of the profile with only three broad Gaussians cannot explain
the observed excess in power at 40 Hz.

Figure 2:  Light curve of the X-ray burst in 18–50 keV obtained by Insight-
HXMT /ME with time resolution of 5 ms. The blue line is the light curve after
dead time correction. The green line represents the fit result of the burst
profile with three Gaussian functions as given in eq. (1). The red line is the
fitted oscillation curve with eq. (11). The two dashed green lines show the
lower and upper envelopes of the amplitude variation of the QPOs described by
eq. (12). The two vertical dashed lines denote the times of the two radio
pulses of FRB 200428. The inset is the light curve around the burst peak with
the solid green line subtracted.  Figure 3:  The Leahy-normalized periodograms
of the light curves of Insight-HXMT/ME with time resolution of 5 ms. The solid
black line represents the Leahy-normalized periodograms. 10000 simulated light
curves sampled from the fitted three-Gaussian burst profile by adding the
Poisson noise is used to generate the Leahy-normalized periodograms in the
same way as we did for the data. Thus the distribution of Leahy power in each
frequency bin can be derived. The mean simulated power in each frequency bin
is shown in red solid line with the maximum (green dots) and minimum (magenta
dots). The blue dotted line means the Leahy power with 2.

The QPO candidate in the light curve as shown in Figure 2 can also be fitted
with,

$Q=N_{\rm g}G(t,t_{\rm g},\sigma_{\rm g})\sin(2\pi{f_{\rm q}(t-t_{\rm
q})})+R_{\rm 3},$ (11)

where $N_{\rm g}$, $t_{\rm g}$ and $\sigma_{\rm g}$ are the normalization,
arrival time and width of the broad Gaussian profile of the X-ray burst.
$f_{\rm q}$ is QPO frequency and $t_{\rm q}$ is the time offset of sinusoidal
function. The fitting results are $N_{\rm g}=1.4\pm 0.2$, $t_{\rm g}=0.452\pm
0.004$, $\sigma_{\rm g}=0.024\pm 0.004$, $f_{\rm q}=40.4\pm 1.2$, $t_{\rm
q}=0.09\pm 0.01$. The frequency $f_{\rm q}$ is consistent with the result from
the periodogram analysis presented below. The lower and upper envelopes of the
amplitude modulation of the QPO candidate can thus be expressed by

$Q_{\rm en}=R_{\rm 3}\pm{N_{\rm g}G(t,t_{\rm g},\sigma_{\rm g})}.$ (12)

The fitted ME light curve with eq. (11) is shown in Figure 2 (b), where the
lower and upper envelopes of the amplitude modulation of the QPO described by
eq. (12) are also plotted.

Table 1: Fitting parameters of the light curve of Insight-HXMT/ME (18-50 keV) with eq.(1). Inst. | $N_{\rm p1}$ | $t_{\rm p1}$ | $\sigma_{\rm p1}$ | $N_{\rm p2}$ | $t_{\rm p2}$ | $\sigma_{\rm p2}$ | $N_{\rm p3}$ | $t_{\rm p3}$ | $\sigma_{\rm p3}$ | $l$ |
---|---|---|---|---|---|---|---|---|---|---|---
ME | 49.7 | 0.446 | 0.044 | 10.4 | 0.545 | 0.095 | 8.02 | 0.246 | 0.065 | 1.129 |

###  3.2 Bayesian method

The Monte Carlo method outlined above is versatile and powerful, but it has
limitations. In this section, we will explore the Bayesian method which is a
conservative method based on the assumption that red noise component dominates
the periodogram. This approach is to model the observed periodogram rather
than the light curve. The data between 0.1 s to 0.74 s is used to generate the
dead time corrected light curve with time resolution of 2.5 ms and the
corresponding periodogram. The 0.64 s duration basically covers the entire
X-ray burst of FRB 200428 as shown in Figure 2. In order to study the non-
stationary effect on this burst, we also compare the results of different
duration time such as the main peak of the burst (0.32 s) and more background
data besides the burst such as 1.28 s.

####  3.2.1 Selection of a noise model

In the following, we briefly illustrate the analysis procedure with the X-ray
burst associated with FRB 200428. This burst has a duration of $T_{90}$ with
0.65 s and the periodogram of 0.64 s (from 0.1 to 0.74 s) length with time
resolution of 2.5 ms for this burst is presented in Figure 4. After fitting
both a simple power law ($H_{0}$ model) and the broken power law ($H_{1}$
model), the likelihood ratio between the two models using eq.(8) is calculated
as LRT = 5.36. The fits to the periodogram and the ratios (data/model) are
presented in Figure 4 using the MAP estimate of the parameters.

Figure 4:  The Leahy-normalized periodogram and residuals of the different
models. Upper panel: The Leahy-normalized periodogram (black), power-law fit
(blue), broken power law fit (red) and broken power law plus the Lorentz fit
(green) of the periodogram. The three lower panels show the smoothed ratio of
the observed periodogram and power law fit, broken power law fit, and broken
power law plus Lorentz fit with MAP estimates, respectively. The smooth is
derived by a Wiener filter with a width of three frequency bins in order to
reduce the probability of the minimization algorithm terminating in local
minima due to sharp noise features. The QPO candidate with center frequency of
40 Hz is clearly shown in the ratio plot.

We utilize the Markov Chain Monte Carlo simulations (MCMCs) implemented in
Python package emcee (Foreman-Mackey et al., 2013) to generate sets of
parameters of posterior probabilities. We then sample from the posterior
distribution of the power law model via MCMC and simulate the periodograms
according to eq. (2) from chains of posterior distribution as shown in Figure
5. These simulated periodograms are again fitted with models $H_{0}$ and
$H_{1}$, respectively, so that we can build the likelihood ratio distribution
of the power law model. Comparing the distribution of simulated likelihood
ratios with that of the observed LRT, we can compute the tail area probability
or the posterior p-value as 0.03 which is shown in Figure 6. Thus we can
reject the power law model and accept the broken power law model for the noise
model.

Figure 5:  The posterior distributions of parameters in the simple power law
are presented on the diagonal. The parameter $A$ is the amplitude and $\alpha$
is the index of the power law as defined in eq.(3). $C$ is the constant
representing the white noise component and also defined in eq.(3). The off-
diagonal panels show correlations between parameters. The scatter plots are
selected from the last 10000 picked parameter pairs from the entire sample of
100000 parameter sets. We can observe that a very tight correlation between
power law index and normalization, and very little correlation between the
normalization and the white noise.  Figure 6:  Distribution of likelihood
ratios for 150000 simulations of the null hypothesis (power-law model). The
observed value of LRT is indicated as the red vertical line. The ratio of the
tail area is about 0.03, indicating that a more complex model (in this case
the broken power law) may be more appropriate in modeling the broadband noise.

####  3.2.2 Searching for QPO

The noise model has been selected above as broken power law, then we can
search for the QPO as a model selection problem. We compare the selected noise
model to a more complex model combining both the noise model and a Lorentz
function to account for the QPO. The burst periodogram are fitted by
bknpow+constant as a simple model ($H_{0}$), and can be also fitted by
bknpow+Lorentz+constant as a nested model ($H_{1}$). bknpow and Lorentz are
the power law function with a break frequency and the Lorentz function in
Xspec. After fitting the periodogram with the two models above, the likelihood
ratio between the two model can be obtained, denoted as
$T_{\rm{LRT}}^{\rm{obs}}$. The $D_{\rm{min}}(H_{0})$ and $D_{\rm{min}}(H_{1})$
of the observed burst periodogram are 743.75 and 712.76, respectively. The
ratio of the periodogram and $H_{1}$ model is also displayed in the bottom
panel of Figure 4. Thus the observed LRT of the burst is
$T_{\rm{LRT}}^{\rm{obs}}=30.99$.

Next, we will clarify how we calibrate the posterior distribution of LRT and
calculate the p-value. Using the $H_{0}$ model, we pick $n$ sets of parameters
for the $H_{0}$ model from the posterior MCMC chains. Each set of parameters
is simulated to produce a periodogram, for which an LRT value for $H_{0}$ and
$H_{1}$ is computed. We then have $n$ LRT values for $n$ periodogram. The
distribution of the $n$ LRT values is the posterior distribution of LRT which
is shown in Figure 8 with red solid line. The observed LRT value with red
vertical line is also shown in Figure 8. As a result, the p-value or the ratio
of the tail areas of $T_{\rm{LRT}}^{\rm{obs}}$ is 1e-6, which corresponds to
4.75 $\sigma$. We find a posterior centroid frequency for the QPO of
$39.3_{-1.20}^{+1.22}$ Hz, with a width of $3.7_{-2.01}^{+3.74}$ Hz. We use
the posterior distribution of QPO parameters to derive the fractional root-
mean-squared (RMS) amplitude $A_{frac}$, and find $A_{frac}=0.36\pm{0.096}$.
The results of 0.64 s duration is also displayed in Table 2.

###  3.3 Significance at different duration

The analysis presented above makes a strong assumption about the data: our
choice of a $\chi^{2}_{2}$ distributed likelihood around the model power
spectrum implies that the periodogram is the result of a stationary noise
process. This assumption is not strictly the case for burst signals, as shown
in Huppenkothen et al. (2013), it is a conservative assumption that holds for
all but the lowest frequencies in the periodogram. In this paper, we also
explore the effects of the non-stationary nature of a transient light curve,
and a QPO within it. We show that selecting more background signals before or
after the burst signal leads to greatly overestimating the significance of
potential QPOs. We also show that the significance will be overestimated if
the whole burst signal is selected for analysis but the QPOs occur in only
part of the burst.

To search for the duration of the 40 Hz QPO, we utilize wavelet analysis and
two most widely used methods: the generalized Lomb-Scargle Periodogram (LSP;
Lomb (1976); Scargle (1982)) and Weighted Wavelet Z-transform (WWZ; Foster
(1996)) are used to obtain the power spectra of the 0.64 s light curve with
time resolution of 2.5 ms. In this work the power spectra of LSP method is
checked with the independent results of WWZ approach. As shown in Figure 7,
the QPO signal appears in the main peak of time range from 0.3 s to 0.6 s.
Then we choose the duration of the main peak about 0.32 s and more background
before and after the main peak, such as 1.28 s to compare their differences.

Using the same method for the 0.64 s, we select the noise model for 0.32 s and
1.28 s. $H_{0}$ model is power law for 0.32 s and broken power law for 1.28 s.
As shown in Table 2, for each duration time, after selecting the noise model
$H_{0}$, we use $H_{0}$ model to generate a large sets of parameters (for
power law, three parameters as shown in eq.(3) and for broken power law, five
parameters as shown in eq.(4)) from their posterior distributions. Each
parameter set is used to simulate a periodogram using the exponential
distribution in eq.(2). For the simulated periodogram we also use $H_{0}$
model and $H_{1}$ model to get the $D_{\rm{min}}$ values as defined in eq.(6)
using MCMC. The $T_{\rm{LRT}}$ can be derived from eq.(8). Then the
distribution of the simulated likelihood ratios can be derived from multiple
sets of parameters of $H_{0}$ model. Comparing the distribution of simulated
likelihood ratios with that of the observed $T_{\rm{LRT}}$ as shown in Figure
8, we can compute the tail area probability and this is also the normalized
posterior p-value. 150000 and 2000000 sets of parameters for the $H_{0}$ model
from the posterior Markov chain are used to simulate the periodogram for 0.32
s and 1.28 s, respectively. The likelihood ratio distributions and the
observed LRT values are also displayed in Figure 8. The posterior parameters
of the QPO, such as the centroid frequency, width and normalization are also
shown in Table 2.

The smaller p-value shown in Table 2 indicates a higher significance of the
QPO. Although the parameters of the QPO are consistent within the
uncertainties in the three durations, however, the uncertainties of the QPO
parameters are smaller for the longer duration. Our explanation about why the
significance is higher for a long duration used is as follows. We find that
the QPO is only present in the main peak of the X-ray burst. However, the
noise level varies through the time series, i.e., the noise level is higher
during the main peak where the count rate is higher (and also varies) than
that outside the main peak. This is a typical non-stationary process.
Extending the selected duration beyond the main peak is equivalent to using
the noise level beyond the main peak to estimate the noise level during the
main peak, which does not change the strength of the QPO but decreases the
estimated noise. Consequently, the signal to noise ratio and detection
significance the QPO are artificially increased.

To test this explanation, we have also done simulations for 0.64 s and 1.28 s
and assume the light curve before and after the main peak (0.32 s) are only
Poisson processes, i.e., by replacing the observed light curve with the
Poisson sampled light curve of the same count rate. We have obtained the same
results: a longer duration outside the main peak leads to a higher
significance of the QPO. Therefore, the reliable way to estimate the
significance of the QPO is to use the segment of data in which the QPO is
present.

As a consequence, we report the final results of QPO from 0.32 s which are
shown in Table 2.

Figure 7:  Wavelet analysis of the ME light curve from 0.1 s to 0.74 s with time resolution of 2.5 ms. The top panel is the light curve and the left bottom panel is the result of the wavelet analysis, the 40 Hz QPO occurs in the time range of 0.3 s to 0.6 s, which is indicated with red dotted lines. The right bottom panel represents the WWZ power and LSP power.  Figure 8:  Distributions of likelihood ratio test for 150000, 2000000 and 2000000 simulations of the null hypothesis for 0.32 s (black), 0.64 s (red) and 1.28 s (blue), respectively. The corresponding $H_{0}$ model is power law, broken power law and broken power law, respectively. The observed values of LRT are indicated as vertical lines. The ratios of the tail area for 0.32 s and 0.64 s are 2.9e-4 and 1e-6, both are indicating that a more complex model including additional Lorentz model may be more appropriate in modeling the periodogram. As for the 1.28 s duration, the tail area is zero and the p-value is an upper limit.  Table 2: Results for the significance of QPO and QPO parameters at different duration time using the LRT statistic. Duration time |  $H_{0}$ model |  $H_{1}$ model | $D_{min}(H_{0})$ | $D_{min}(H_{1})$ | $T_{\rm{LRT}}^{\rm{obs}}$ | p-value | Frequency | Width | Norm
---|---|---|---|---|---|---|---|---|---
(s) |  |  |  |  |  |  | (Hz) | (Hz) |
0.32 | PL | PL+Lorentz | 378.91 | 359.07 | 19.84 | 2.9e-4 | $39.14_{-1.94}^{+1.88}$ | $4.88_{-2.73}^{+5.67}$ | $389.15_{-167.17}^{+260.55}$
0.64 | BPL | BPL+Lorentz | 743.75 | 712.76 | 30.99 | 1e-6 | $39.28_{-1.20}^{+1.22}$ | $3.67_{-2.01}^{+3.74}$ | $315.27_{-116.74}^{+218.63}$
1.28 | BPL | BPL+Lorentz | 1479.71 | 1421.63 | 58.08 | <5e-7 | $39.26_{-0.8}^{+0.78}$ | $3.24_{-1.41}^{+2.02}$ | $267.95_{-76.85}^{+135.87}$

As described in the section 2.3.2, another statistics for investigating narrow
QPOs is $T_{\rm{R}}$, which is defined in eq.(9). We also search the
periodogram for the highest data/model outlier and compare this outlier to
those distributed by pure noise to find narrow features that may be candidates
for a possible QPO. One shortcoming of the $T_{\rm{R}}$ statistic is that it
optimally detects periodic signals confined to one frequency bin. Since the
QPO power spreads over several bins, this is not an optimal way of detecting
broad signals. There are several ways to reduce this restriction. One is to
bin (or smooth) the data in some way, and compute $T_{\rm{R}}$ for the binned
data. If we bin the simulated periodograms in the same way, then the test
statistic $T_{\rm{R}}$ for the binned data is comparable to the distribution
approximated by our simulations, and the latter can be used to derive
posterior predictive p-values.

We use the same chosen noise model ($H_{0}$ model) to fit the periodogram and
same samples of simulated periodograms in the previous step to search for the
highest data/model outlier. We bin the periodogram geometrically, where the
bin size grows with frequency and the bin factor is shown in Table 3. Figure 9
also presents as an example for the bin of 1.28 s periodogram and we choose
the bin factor as 0.25. We also compare the observed value of the $T_{\rm{R}}$
statistic and compute the corresponding p-values for the three duration of the
burst. The results are summarized in Table 3. The QPO with central frequency
at about 40 Hz can be found for different duration. The logarithmic rebin of
the periodogram can not ensure all the QPO frequency bins are in the same bin,
therefore, this method are just used as a verification of the QPO.

Figure 9:  The normalized Leahy power of ME with the light curve of 1.28 s. The black line is the raw PDS of ME, the blue and red lines are the logarithmic rebining of the raw PDS with factor of 0.2 and 0.25, respectively. The rebin factor of 0.25 is more suitable for the QPO signal.  Table 3: Results for the significance of QPO and QPO parameters at different duration time using the $T_{\rm{R}}$ statistic. Duration time(s) |  $H_{0}$ model | rebin factor | Frequency (Hz) | $T_{\rm{R}}$ | p-value
---|---|---|---|---|---
0.32 | PL | 0.2 | 37.3 | 10.6 | 0.0047
0.64 | BPL | 0.2 | 37.3 | 10.3 | 0.00096
1.28 | BPL | 0.25 | 38.2 | 9.8 | 0.00024

###  3.4 Cross correlation analysis

To further verify the robustness of the $\sim$40 Hz QPO detection, cross-
correlations between the detrended light curves of Insight-HXMT/ME and Konus-
Wind are made and shown in Figure 10, along with cross-correlations between
each light curve and white noise. The QPOs in both light curves match each
other both in phase and period of 24.5 ms, while no statistically meaningful
correlation is found with white noise.

Figure 10:  Cross correlation between the light curves of Insight-HXMT/ME
(18–50 keV) and KW (18–80 keV). Blue points and solid line: coefficient of the
cross correlation between the detrended light curves of Insight-HXMT/ME and
KW. Green points and line: between the detrended light curve of Insight-
HXMT/ME and white noise. Red points and line: between the detrended light
curve of KW and white noise. The errors are estimated with Monte Carlo
simulations. The vertical dashed lines mark the QPO period of $\sim$$24.5$ ms.

Cross correlation function (CCF) is used here to understand if the QPO signals
found in the light curves of Insight-HXMT/ME and Konus-Wind are synchronised.
The CCF between the detrended light curves (after removing the three broad
Gaussian peaks by fitting with eq. (1) of ME and Konus-Wind shows strong
oscillations of $\sim 24.5$ ms from $\tau=-0.1$ s to $\tau=0.05$ s as plotted
in Figure 10. The errors for CCF are estimated from Monte-Carlo method, by
sampling the light curves 1000 times. We also make the cross correlation
between the detrended light curves and white noise with the same average
counting rate, as expected, there are no oscillation structures in CCFs.
However, there are some oscillation signals around $-0.15$ s and $0.15$ s,
which are also generated by the oscillation structures in the real light
curves.

To further understand the above results of cross correlation analysis, we use
sinusoidal signals to simulate the structures in CCFs obtained above. As shown
in Figure 11 (upper panel: background + sinusoidal signals), two sinusoidal
signals of amplitudes of 8 and 24 are simulated with background 24, which have
frequency $10/2\pi$ Hz and only appear for 3 seconds in the light curves. The
Poisson fluctuations are added on the light curves. The background light
curves only including white noise are generated from Poisson sampling with
expected count of 24 in each time bin, as shown in Figure 11 (lower panel:
background only) for two of the simulated background light curves of white
noise.

We first make cross correlations between the two simulated light curves
(including their backgrounds) with different amplitudes of sinusoidal
oscillations. Very significant and periodic structures are shown in Figure 12,
demonstrating the power of cross correlation analysis in identifying weakly
oscillating signals (‘Seq 2’: blue line in the upper panel of Figure 11),
given a template of the oscillations (‘Seq 1’: red line in the upper panel of
Figure 11). Therefore, the cross correlation analysis presented in Figure 10
is equivalent to using ME’s light curve, which contains strong QPO signals, as
a template in identifying the relatively weak QPOs in the light curve of
Konus-Wind.

The CCFs between the white noise and both light curves with weak and strong
oscillation signals are shown in Figure 13. The upper panel shows the result
from one single simulation run, where many oscillatory structures exist
(Figure 13 upper panel). The processes are repeated for 1000 times to obtain
the mean CCFs, which show clear oscillation structures around $\pm 5$ (lower
panel in Figure 13). These structures are similar and should have the same
origin as the structures near both ends of the correlation functions of white
noise with the light curves of ME and Konus-Wind, as shown in Figure 10.

Figure 11:  Simulated light curves with sinusoidal signals and only white
noise. Upper panel (a): simulated light curves with 3 seconds of sinusoidal
signals of amplitudes 8 and 24, respectively. Lower panel (b): simulated light
curves with only white noise with Poisson expectation value of 24 in each time
bin.  Figure 12:  Cross correlations between the two simulated light curves
with strong and weak sinusoidal oscillations. Upper panel (a): One cross
correlation for the two simulated light curves (blue and red lines in the
upper panel of Figure 11). Lower panel (b): mean cross correlations for the
1000 simulated light curves. Therefore the strong sinusoidal oscillations are
used as a template to identify effectively the weak sinusoidal oscillations,
which can only be seen marginally alone.  Figure 13:  The cross correlations
for the simulated light curves without sinusoidal oscillations. Upper panel
(a): One cross correlations for the simulated light curves. Lower panel (b):
mean cross correlations for the 1000 simulated light curves.

##  4 Physical interpretation

###  4.1 Constraints on theoretical models

The 40 Hz QPO reported in this paper poses the following constraints on the
available models for FRB 200428:

  * •

The fact that this burst is special with respect to the majority of SGR bursts
suggest that the rarity of FRB-SGR burst associations (Lin et al., 2020d) is
likely due to intrinsic rather than extrinsic (beaming or narrow spectra; Lin
et al. (2020d)) factors.

  * •

Since QPOs are directly related with crust oscillations, and since the two
spikes on the X-ray lightcurve that are associated with the two pulses of FRB
200428 are roughly consistent with QPO peaks, it is almost certain that the
X-ray emission as well as the associated FRB pulses originate from the
magnetosphere of the neutron star (Lu et al., 2020; Wang, 2020; Wang et al.,
2020; Wadiasingh & Chirenti, 2020; Katz, 2020; Yang et al., 2020). Models
invoking emission outside the magentosphere (Margalit et al., 2020; Yu et al.,
2020) are thus disfavored, unless the ejections of the emitter are associated
with crust activities (Yuan et al., 2020).

  * •

The trigger of QPOs in the neutron star crust is likely from internal, but an
external trigger powered by infall of an asteroid is also possible (Dai, 2020;
Geng et al., 2020), which might explain the rarity of the FRB-associated
bursts. Such models, however, need to interpret the existence of QPOs already
before the spikes associated with FRB pulses.

###  4.2 A plausible QPO model

The identified QPO at 40 Hz may be explained as the toroidal mode of the crust
(Duncan, 1998), or the magnetoelastic modes if the coupling effect between the
crust and core is considered (van Hoven & Levin, 2011). For the crustal
toroidal mode, the 40 Hz QPO can be explained as the ${}_{3}t_{0}$ mode if the
tangled magnetic field in the crust is $B_{t}\ll B_{\mu}$, or the fundamental
toroidal mode ${}_{2}t_{0}$ if $B_{t}\approx B_{\mu}$, where ${}_{\ell}t_{n}$
is the notation for the general toroidal modes, $\ell$ is the angular quantum
number, and $n$ is the number of radial nodes in the eigenfunctions, and
$B_{\mu}=(4\pi\mu)^{1/2}\approx 4\times 10^{15}\rho_{14}^{0.4}$ G with $\mu$
being the crust shear modulus, and $\rho_{14}$ being crust density in units of
$10^{14}$g cm-3 (Duncan, 1998). Note the internal magnetic field is usually
expected to be much higher than the surface dipole magnetic field, which is
$2.2\times 10^{14}$ G for SGR J1935+2154 (Israel et al., 2016).

##  5 Discussion and Summary

Magnetar bursts are a potential window for the interior of neutron stars, via
the oscillations measured in magnetar giant flares. As for the X-ray burst
associated with FRB 200428, different strategies are used to search for the
QPO and confirm the significance of the QPO. Monte Carlo simulations of light
curves fail to be predictive when there is no precise knowledge of the burst
profile. Moreover, there is a degeneracy between the burst profile, a
potential red noise component, and the QPO.

In the absence of the burst profile, we advocate a conservative Bayesian
method which models the Leahy-normalized power of the burst as a pure red
noise process. The conventional method to generate Leahy-normalized power is
used the whole light curves of the X-ray burst. However, when the QPO is
present in part of the burst phase, we may overestimate the significance of
the QPO signal and underestimate the errors of QPO parameters. The reliable
way to estimate the significance of the QPO and its parameters in non-
stationary processes is to use the time segment of data in which the QPO is
present.

In the peculiar X-ray burst associated with FRB 200428 from the Galactic
magnetar SGR J1935+2154, we have discovered a QPO with a central frequency of
$39.14_{-1.94}^{+1.88}$ Hz and width of $4.88_{-2.73}^{+5.67}$ between 18–50
keV with p-value of 2.9e-4. The frequency, phase, modulation pattern and
energy range all point to a physical connection with the two X-ray spikes
coinciding with the two radio pulses of FRB 200428. The rarity of both QPOs in
magnetar X-ray bursts and short radio pulses from magnetars indicates that the
rarity of FRB-SGR burst associations (Lin et al., 2020d) is very likely due to
intrinsic reasons. It also poses important constraints on available FRB models
proposed to interpret FRB 200428. In particular, models invoking relativistic
shocks (Margalit et al., 2020) are disfavored. Instead, our QPO result
supports these models invoking crust oscillation as the driver of
magnetospheric emission from a magnetar magnetosphere (Lu et al., 2020; Wang,
2020; Wadiasingh & Timokhin, 2019). The specific frequency (40 Hz) of the QPOs
may also pose interesting constraints on the magnetic field configurations in
the neutron star crust.

##  6 Acknowledgments

Acknowledgements This work made use of the data from the Insight-HXMT mission,
a project funded by China National Space Administration (CNSA) and the Chinese
Academy of Sciences (CAS). We gratefully acknowledge the support from the
National Program on Key Research and Development Project (Grant
No.2021YFA0718500) from the Minister of Science and Technology of China
(MOST). The authors thank supports from the National Natural Science
Foundation of China under Grants U1838201, U1838202, U1938109, U1938102,
U1938108, 11473027, 11733009, 1173309, Y829113, 11673023, 11703002, 11833003.
We also thank the support from the Strategic Priority Program on Space
Science, the Chinese Academy of Sciences, Grant No.XDA15020503. J.S.Wang
acknowledges the support from the Alexander von Humboldt Foundation.

## References

  * Bochenek et al. (2020) Bochenek, C. D., Ravi, V., Belov, K. V., et al. 2020, Nature, 587, 59, doi: 10.1038/s41586-020-2872-x
  * Cao et al. (2020) Cao, X., Jiang, W., Meng, B., et al. 2020, Science China Physics, Mechanics, and Astronomy, 63, 249504, doi: 10.1007/s11433-019-1506-1
  * Castro-Tirado et al. (2021) Castro-Tirado, A. J., Østgaard, N., Göǧüş, E., et al. 2021, Nature, 600, 621, doi: 10.1038/s41586-021-04101-1
  * Chen et al. (2020) Chen, Y., Cui, W., Li, W., et al. 2020, Science China Physics, Mechanics, and Astronomy, 63, 249505, doi: 10.1007/s11433-019-1469-5
  * CHIME/FRB Collaboration et al. (2020) CHIME/FRB Collaboration, Andersen, B. C., Bandura, K. M., et al. 2020, Nature, 587, 54, doi: 10.1038/s41586-020-2863-y
  * Dai (2020) Dai, Z. G. 2020, arXiv e-prints, arXiv:2005.12048.  https://arxiv.org/abs/2005.12048
  * Duncan (1998) Duncan, R. C. 1998, ApJ, 498, L45, doi: 10.1086/311303
  * Foreman-Mackey et al. (2013) Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, Publications of the Astronomical Society of the Pacific, 125, 306
  * Foster (1996) Foster, G. 1996, AJ, 111, 541, doi: 10.1086/117805
  * Fox et al. (2001) Fox, D. W., Lewin, W. H. G., Rutledge, R. E., et al. 2001, MNRAS, 321, 776, doi: 10.1046/j.1365-8711.2001.04085.x
  * Geng et al. (2020) Geng, J.-J., Li, B., Li, L.-B., et al. 2020, arXiv e-prints, arXiv:2006.04601.  https://arxiv.org/abs/2006.04601
  * Guidorzi (2011) Guidorzi, C. 2011, MNRAS, 415, 3561, doi: 10.1111/j.1365-2966.2011.18969.x
  * Huppenkothen et al. (2014a) Huppenkothen, D., Heil, L. M., Watts, A. L., & Göğüş, E. 2014a, ApJ, 795, 114, doi: 10.1088/0004-637X/795/2/114
  * Huppenkothen et al. (2014b) Huppenkothen, D., Watts, A. L., & Levin, Y. 2014b, ApJ, 793, 129, doi: 10.1088/0004-637X/793/2/129
  * Huppenkothen et al. (2013) Huppenkothen, D., Watts, A. L., Uttley, P., et al. 2013, ApJ, 768, 87, doi: 10.1088/0004-637X/768/1/87
  * Huppenkothen et al. (2014c) Huppenkothen, D., D’Angelo, C., Watts, A. L., et al. 2014c, ApJ, 787, 128, doi: 10.1088/0004-637X/787/2/128
  * Huppenkothen et al. (2019) Huppenkothen, D., Bachetti, M., Stevens, A. L., et al. 2019, ApJ, 881, 39, doi: 10.3847/1538-4357/ab258d
  * Israel et al. (2005) Israel, G. L., Belloni, T., Stella, L., et al. 2005, ApJ, 628, L53, doi: 10.1086/432615
  * Israel et al. (2016) Israel, G. L., Esposito, P., Rea, N., et al. 2016, MNRAS, 457, 3448, doi: 10.1093/mnras/stw008
  * Katz (2020) Katz, J. I. 2020, arXiv e-prints, arXiv:2006.03468.  https://arxiv.org/abs/2006.03468
  * Kirsten et al. (2021) Kirsten, F., Snelders, M. P., Jenkins, M., et al. 2021, Nature Astronomy, 5, 414, doi: 10.1038/s41550-020-01246-3
  * Kouveliotou et al. (1998) Kouveliotou, C., Dieters, S., Strohmayer, T., et al. 1998, Nature, 393, 235, doi: 10.1038/30410
  * Li et al. (2021) Li, C. K., Lin, L., Xiong, S. L., et al. 2021, Nature Astronomy, 5, 378\.  https://arxiv.org/abs/2005.11071
  * Li et al. (2020) Li, X., Li, X., Tan, Y., et al. 2020, Journal of High Energy Astrophysics, 27, 64, doi: 10.1016/j.jheap.2020.02.009
  * Lin et al. (2020a) Lin, L., Göğüş, E., Roberts, O. J., et al. 2020a, ApJ, 902, L43, doi: 10.3847/2041-8213/abbefe
  * Lin et al. (2020b) Lin, L., Göğüs, E., Roberts, O. J., et al. 2020b, ApJ, 893, 156, doi: 10.3847/1538-4357/ab818f
  * Lin et al. (2020c) Lin, L., Zhang, C. F., Wang, P., et al. 2020c, Nature, 587, 63, doi: 10.1038/s41586-020-2839-y
  * Lin et al. (2020d) —. 2020d, arXiv e-prints, arXiv:2005.11479.  https://arxiv.org/abs/2005.11479
  * Liu et al. (2020) Liu, C., Zhang, Y., Li, X., et al. 2020, Science China Physics, Mechanics, and Astronomy, 63, 249503, doi: 10.1007/s11433-019-1486-x
  * Lomb (1976) Lomb, N. R. 1976, Ap&SS, 39, 447, doi: 10.1007/BF00648343
  * Lorimer et al. (2007) Lorimer, D. R., Bailes, M., McLaughlin, M. A., Narkevic, D. J., & Crawford, F. 2007, Science, 318, 777, doi: 10.1126/science.1147532
  * Lu et al. (2020) Lu, W., Kumar, P., & Zhang, B. 2020, arXiv e-prints, arXiv:2005.06736.  https://arxiv.org/abs/2005.06736
  * Margalit et al. (2020) Margalit, B., Beniamini, P., Sridhar, N., & Metzger, B. D. 2020, arXiv e-prints, arXiv:2005.05283.  https://arxiv.org/abs/2005.05283
  * Mereghetti et al. (2020) Mereghetti, S., Savchenko, V., Ferrigno, C., et al. 2020, ApJ, 898, L29, doi: 10.3847/2041-8213/aba2cf
  * Miller et al. (2019) Miller, M. C., Chirenti, C., & Strohmayer, T. E. 2019, ApJ, 871, 95, doi: 10.3847/1538-4357/aaf5ce
  * Protassov et al. (2002) Protassov, R., van Dyk, D. A., Connors, A., Kashyap, V. L., & Siemiginowska, A. 2002, ApJ, 571, 545, doi: 10.1086/339856
  * Ridnaia et al. (2021) Ridnaia, A., Svinkin, D., Frederiks, D., et al. 2021, Nature Astronomy, 5, 372, doi: 10.1038/s41550-020-01265-0
  * Scargle (1982) Scargle, J. D. 1982, ApJ, 263, 835, doi: 10.1086/160554
  * Strohmayer & Watts (2005) Strohmayer, T. E., & Watts, A. L. 2005, ApJ, 632, L111, doi: 10.1086/497911
  * Strohmayer & Watts (2006) —. 2006, ApJ, 653, 593, doi: 10.1086/508703
  * Tavani et al. (2021) Tavani, M., Casentini, C., Ursi, A., et al. 2021, Nature Astronomy, 5, 401, doi: 10.1038/s41550-020-01276-x
  * Thompson & Duncan (1995) Thompson, C., & Duncan, R. C. 1995, MNRAS, 275, 255, doi: 10.1093/mnras/275.2.255
  * Timmer & Koenig (1995) Timmer, J., & Koenig, M. 1995, A&A, 300, 707
  * Tuo et al. (2021) Tuo, Y., Li, X., Ge, M., et al. 2021, arXiv e-prints, arXiv:2109.04709.  https://arxiv.org/abs/2109.04709
  * Van der Klis (1989) Van der Klis, M. 1989, in Timing neutron stars (Springer), 27–69
  * van Hoven & Levin (2011) van Hoven, M., & Levin, Y. 2011, MNRAS, 410, 1036, doi: 10.1111/j.1365-2966.2010.17499.x
  * Vaughan (2010) Vaughan, S. 2010, MNRAS, 402, 307, doi: 10.1111/j.1365-2966.2009.15868.x
  * Wadiasingh & Chirenti (2020) Wadiasingh, Z., & Chirenti, C. 2020, arXiv e-prints, arXiv:2006.16231.  https://arxiv.org/abs/2006.16231
  * Wadiasingh & Timokhin (2019) Wadiasingh, Z., & Timokhin, A. 2019, ApJ, 879, 4, doi: 10.3847/1538-4357/ab2240
  * Wang (2020) Wang, J.-S. 2020, arXiv e-prints, arXiv:2006.14503.  https://arxiv.org/abs/2006.14503
  * Wang et al. (2021) Wang, P., Cai, C., Zhao, X. Y., et al. 2021, GRB Coordinates Network, 29377, 1
  * Wang et al. (2020) Wang, W.-Y., Xu, R., & Chen, X. 2020, arXiv e-prints, arXiv:2005.02100.  https://arxiv.org/abs/2005.02100
  * Xiao et al. (2022) Xiao, S., Zhao, Y., Zhang, Y. Q., et al. 2022, GRB Coordinates Network, 31443, 1
  * Yang et al. (2020) Yang, Y.-P., Zhu, J.-P., Zhang, B., & Wu, X.-F. 2020, arXiv e-prints, arXiv:2006.03270.  https://arxiv.org/abs/2006.03270
  * Younes et al. (2017) Younes, G., Kouveliotou, C., Jaodand, A., et al. 2017, ApJ, 847, 85, doi: 10.3847/1538-4357/aa899a
  * Younes et al. (2021) Younes, G., Baring, M. G., Kouveliotou, C., et al. 2021, Nature Astronomy, 5, 408.  https://arxiv.org/abs/2006.11358
  * Yu et al. (2020) Yu, Y.-W., Zou, Y.-C., Dai, Z.-G., & Yu, W.-F. 2020, arXiv e-prints, arXiv:2006.00484.  https://arxiv.org/abs/2006.00484
  * Yuan et al. (2020) Yuan, Y., Beloborodov, A. M., Chen, A. Y., & Levin, Y. 2020, arXiv e-prints, arXiv:2006.04649.  https://arxiv.org/abs/2006.04649
  * Zhang et al. (2020a) Zhang, C. F., Jiang, J. C., Men, Y. P., et al. 2020a, The Astronomer’s Telegram, 13699, 1
  * Zhang et al. (2020b) Zhang, S. N., Xiong, S. L., Li, C. K., et al. 2020b, The Astronomer’s Telegram, 13696, 1
  * Zhang et al. (2020c) Zhang, S.-N., Li, T., Lu, F., et al. 2020c, Science China Physics, Mechanics, and Astronomy, 63, 249502, doi: 10.1007/s11433-019-1432-6