text:
\section{Introduction}
Neural Machine Translation (NMT) has advanced the state of the art in MT~\citep{NIPS2014_5346, DBLP:journals/corr/BahdanauCB14, NIPS2017_7181}, but is susceptible to domain shift.
\citet{koehn-knowles-2017-six} consider out-of-domain translation one of the key challenges in NMT.
Such translations may be fluent, but completely unrelated to the input (\textit{hallucinations}), and their misleading nature makes them particularly problematic.

We hypothesise that \textit{exposure bias}~\citep{DBLP:journals/corr/RanzatoCAZ15}, a discrepancy between training and inference, makes this problem worse.
Specifically, training with teacher forcing only exposes the model to gold history, while previous predictions during inference may be erroneous.
Thus, the model trained with teacher forcing may over-rely on previously predicted words, which would exacerbate error propagation.
Previous work has sought to reduce exposure bias in training~\citep{NIPS2015_5956, DBLP:journals/corr/RanzatoCAZ15, shen-etal-2016-minimum, wiseman-rush-2016-sequence, zhang-etal-2019-bridging}.
However, the relevance of error propagation is under debate: \citet{wu-etal-2018-beyond} argue that its role is overstated in literature, and that linguistic features explain some of the accuracy drop at higher time steps.

Previous work has established a link between domain shift and hallucination in NMT \citep{koehn-knowles-2017-six,mller2019domain}.
In this paper, we will aim to also establish an empirical link between hallucination and exposure bias.
Such a link will deepen our understanding of the hallucination problem, but also has practical relevance, e.g.\ to help predicting in which settings the use of sequence-level objectives is likely to be helpful.
We further empirically confirm the link between exposure bias and the `beam search problem', i.e.\ the fact that translation quality does not increase consistently with beam size~\citep{koehn-knowles-2017-six, ott2018analyzing, stahlberg-byrne-2019-nmt}.

We base our experiments on German$\to$English IWSLT'14, and two datasets used to investigate domain robustness by \citet{mller2019domain}: a selection of corpora from OPUS~\citep{lison-tiedemann-2016-opensubtitles2016} for German$\to$English, and a low-resource German$\to$Romansh scenario.
We experiment with Minimum Risk Training (MRT)~\citep{och-2003-minimum, shen-etal-2016-minimum}, a training objective which inherently avoids exposure bias.

Our experiments show that MRT indeed improves quality more in out-of-domain settings, and reduces the amount of hallucination.
Our analysis of translation uncertainty also shows how the MLE baseline over-estimates the probability of random translations at all but the initial time steps, and how MRT mitigates this problem.
Finally, we show that the beam search problem is reduced by MRT.

\section{Minimum Risk Training}

The de-facto standard training objective in NMT is to minimize the negative log-likelihood $\mathcal{L}(\boldsymbol{\theta})$ of the training data $D$\footnote{This is equivalent to maximizing the likelihood of the data, hence \textit{Maximum Likelihood Estimation} (MLE).}:
\begin{equation}
\label{eq:1}
\mathcal{L}(\boldsymbol{\theta})=\sum_{(\mathbf{x},\mathbf{y}) \in D} \sum_{t=1}^{|\mathbf{y}|} -\log P\left(\mathbf{y}_{t} | \mathbf{x}, \mathbf{y}_{<t} ; \boldsymbol{\theta}\right)
\end{equation}
where $\mathbf{x}$ and $\mathbf{y}$ are the source and target sequence, respectively, $\mathbf{y}_t$ is the t\textsuperscript{th} token in $\mathbf{y}$, and $\mathbf{y}_{<t}$ denotes all previous tokens. MLE is typically performed with teacher forcing, where $\mathbf{y}_{<t}$ are ground-truth labels in training, which creates a mismatch to inference, where $\mathbf{y}_{<t}$ are model predictions.

Minimum Risk Training (MRT) is a sequence-level objective that avoids this problem. Specifically, the objective function of MRT is the expected loss (\textit{risk}) with respect to the posterior distribution:
\begin{equation}
\label{eq:2}
\mathcal{R}(\boldsymbol{\theta})=\sum_{(\mathbf{x},\mathbf{y}) \in D} \sum_{\mathbf{\tilde{y}} \in \mathcal{Y}(\mathbf{x})}P\left(\mathbf{\tilde{y}} | \mathbf{x} ; \boldsymbol{\theta}\right) \Delta\left(\mathbf{\tilde{y}}, \mathbf{y}\right)
\end{equation}
in which the loss $\Delta\left(\mathbf{\tilde{y}}, \mathbf{y}\right)$ indicates the discrepancy between the gold translation $\mathbf{y}$ and the model prediction $\mathbf{\tilde{y}}$. Due to the intractable search space, the posterior distribution $\mathcal{Y}(\mathbf{x})$ is approximated by a subspace $\mathcal{S}(\mathbf{x})$ by sampling a certain number of candidate translations, and normalizing:
\begin{equation}
\label{eq:3}
\tilde{P}\left(\mathbf{\tilde{y}} | \mathbf{x} ; \boldsymbol{\theta}, \alpha\right)=\frac{P\left(\mathbf{\tilde{y}} | \mathbf{x} ; \boldsymbol{\theta}\right)^{\alpha}}{\sum_{\mathbf{y}^{\prime} \in \mathcal{S}\left(\mathbf{x}\right)} P\left(\mathbf{y}^{\prime} | \mathbf{x} ; \boldsymbol{\theta}\right)^{\alpha}}
\end{equation}
where $\alpha$ is a hyperparameter to control the sharpness of the subspace.
Based on preliminary results, we use random sampling to generate candidate translations, and following~\citet{edunov-etal-2018-classical}, do not add the reference translation to the subspace.

\section{Experiments}
\subsection{Data}
To verify the effectiveness of our MRT implementation on top of a strong Transformer baseline \citep{NIPS2017_7181}, we first conduct experiments on IWSLT'14 German$\to$English (DE$\to$EN) \citep{cettolo2014report}, which consists of \num{180000} sentence pairs. We follow previous work for data splits~\citep{DBLP:journals/corr/RanzatoCAZ15,edunov-etal-2018-classical}.

For experiments with domain shift, we use data sets and preprocessing as~\citet{mller2019domain}\footnote{\url{https://github.com/ZurichNLP/domain-robustness}}.
For DE$\to$EN, data comes from OPUS~\citep{lison-tiedemann-2016-opensubtitles2016}, and is comprised of five domains: \textit{medical}, \textit{IT}, \textit{law}, \textit{koran} and \textit{subtitles}.
We use \textit{medical} for training and development, and report results on an in-domain test set and the four other domains (out-of-domain; OOD).
German$\to$Romansh (DE$\to$RM) is a low-resource language pair where robustness to domain shift is of practical relevance. The training data is from the Allegra corpus~\citep{scherrer-cartoni-2012-trilingual} (\textit{law} domain) with \num{100000} sentence pairs. The test domain are \textit{blogs}, using data from Convivenza\footnote{\url{https://www.suedostschweiz.ch/blogs/convivenza}}.
We have access to 2000 sentences for development and testing, respectively, in each domain.

We tokenise and truecase data sets with Moses~\citep{koehn-etal-2007-moses}, and use shared BPE with \num{32000} units~\citep{sennrich-etal-2016-neural}.

\subsection{Model}
We implement\footnote{Code available at \url{https://github.com/zippotju/Exposure-Bias-Hallucination-Domain-Shift}} MRT in the Nematus toolkit~\citep{sennrich-etal-2017-nematus}.
All our experiments use the Transformer architecture~\citep{NIPS2017_7181}.
Following ~\citet{edunov-etal-2018-classical}, we use \mbox{$1$-$\text{BLEU}_\text{smooth}$} \citep{lin-och-2004-orange} as the MRT loss.
Models are pre-trained with the token-level objective MLE and then fine-tuned with MRT.
Hyperparameters mostly follow previous work \citep{edunov-etal-2018-classical,mller2019domain}; for MRT, we conduct limited hyperparameter search on the IWSLT'14 development set, including learning rate, batch size, and the sharpness parameter $\alpha$.
We set the number of candidate translations for MRT to 4 to balance effectiveness and efficiency.
Detailed hyperparameters are reported in the Appendix.

\subsection{Evaluation}

For comparison to previous work, we report lowercased, tokenised BLEU~\citep{papineni-etal-2002-bleu} with \textit{multi-bleu.perl} for IWSLT'14, and cased, detokenised BLEU with SacreBLEU~\citep{post-2018-call}\footnote{Signature: BLEU+c.mixed+\#.1+s.exp+tok.13a+v.1.4.2} otherwise.
For settings with domain shift, we report average and standard deviation of 3 independent training runs to account for optimizer instability.

The manual evaluation was performed by two native speakers of German who completed bilingual (German/English) high school or University programs.
We collected $\sim$3600 annotations in total, spread over 12 configurations.
We ask annotators to evaluate translations according to fluency and adequacy. For fluency, the annotator classifies a translation as fluent, partially fluent or not fluent; for adequacy, as adequate, partially adequate or inadequate.
We report kappa coefficient ($K$)~\citep{carletta-1996-assessing} for inter-annotator and intra-annotator agreement in Table~\ref{tab:7}, and assess statistical significance with Fisher's exact test (two-tailed).

\begin{table}
\centering
\begin{small}
\setlength{\tabcolsep}{0.5em}
\begin{tabular}{lcccccc}
\toprule
& \multicolumn{3}{c}{\textbf{inter-annotator}} & \multicolumn{3}{c}{\textbf{intra-annotator}} \\
\cmidrule(lr){2-4}\cmidrule(lr){5-7}
annotation & $P(A)$ & $P(E)$ & $K$ & $P(A)$ & $P(E)$ & $K$ \\
\midrule
fluency & 0.66 & 0.38 & 0.44 & 0.87 & 0.42 & 0.77 \\
adequacy & 0.82 & 0.61 & 0.54 & 0.93 & 0.66 & 0.79 \\
\bottomrule
\end{tabular}
\end{small}
\caption{Inter-annotator (N=250) and intra-annotator agreement (N=617) of manual evaluation.}
\label{tab:7}
\end{table}

\subsection{Results}

Table~\ref{tab:1} shows results for IWSLT'14.
We compare to results by \citet{edunov-etal-2018-classical}, who use a convolutional architecture~\citep{Gehring:2017:CSS:3305381.3305510}, and \citet{wu2018pay}, who report results with Transfomer-base and dynamic convolution.
\begin{table}
\centering
\begin{small}
\begin{tabular}{lc}
\toprule
system   & BLEU \\
\midrule
ConvS2S (MLE)~\citep{edunov-etal-2018-classical} & 32.2 \phantom{(+\textbf{0.6})} \\
ConvS2S (MRT)~\citep{edunov-etal-2018-classical} & 32.8 (+\textbf{0.6}) \\
Transformer (MLE)~\citep{wu2018pay} & 34.4 \phantom{(+\textbf{0.6})} \\
DynamicConv (MLE)~\citep{wu2018pay} & 35.2 \phantom{(+\textbf{0.6})} \\
\midrule
MLE  & 34.7 \phantom{(+\textbf{0.6})}  \\
MRT & 35.2 (+\textbf{0.5}) \\
\bottomrule
\end{tabular}
\end{small}
\caption{Results for IWSLT'14 DE$\to$EN with MLE and MRT (in brackets, improvement over MLE).}
\label{tab:1}
\end{table}

\begin{table*}[!t]
\centering
\begin{small}
\begin{tabular}{lcccc}
\toprule
& \multicolumn{2}{c}{\textbf{DE$\to$EN}} & \multicolumn{2}{c}{\textbf{DE$\to$RM}} \\
\cmidrule(lr){2-3}\cmidrule(lr){4-5}
system & in-domain & average OOD & in-domain & average OOD\\
\midrule
SMT~\cite{mller2019domain} & 58.4 \phantom{($\pm$0.00)} & 11.8 \phantom{($\pm$0.00)}& 45.2 \phantom{($\pm$0.00)}& 15.5 \phantom{($\pm$0.00)}\\
NMT~\cite{mller2019domain} & 61.5 \phantom{($\pm$0.00)} & 11.7 \phantom{($\pm$0.00)} & 52.5 \phantom{($\pm$0.00)} & 18.9 \phantom{($\pm$0.00)}\\
NMT+RC+SR+NC~\cite{mller2019domain} & 60.8 \phantom{($\pm$0.00)} & 13.1 \phantom{($\pm$0.00)} & 52.4 \phantom{($\pm$0.00)} & 20.7 \phantom{($\pm$0.00)}\\
\midrule
MLE w/o LS  & 58.3 ($\pm$0.53) & \phantom{0}9.7 ($\pm$0.25)  & 52.2 ($\pm$0.19) & 15.8 ($\pm$0.39)\\
+MRT & 58.4 ($\pm$0.39) & 10.2 ($\pm$0.26) & 52.1 ($\pm$0.08) & 15.9 ($\pm$0.28)\\
MLE w/ LS & 58.9 ($\pm$0.45) & 11.2 ($\pm$0.16) & 53.9 ($\pm$0.16) & 18.0 ($\pm$0.17)\\
+MRT & 58.8 ($\pm$0.36) & 12.0 ($\pm$0.29) & 53.9 ($\pm$0.12) & 18.7 ($\pm$0.09)\\
\bottomrule
\end{tabular}
\end{small}
\caption{Average BLEU and standard deviation on in-domain and out-of-domain test sets for models trained on OPUS (DE$\to$EN) and Allegra (DE$\to$RM). RC: reconstruction; SR: subword regularization, NC: noisy channel.}
\label{tab:2}
\end{table*}

With 34.7 BLEU, our baseline is competitive.
We observe an improvement of 0.5 BLEU from MRT, comparable to \citet{edunov-etal-2018-classical}, although we start from a stronger baseline (+2.5 BLEU).

Table~\ref{tab:2} shows results for data sets with domain shift.
To explore the effect of label smoothing~\citep{7780677}, we train baselines with and without label smoothing.
MLE with label smoothing performs better by itself, and we also found MRT to be more effective on top of the initial model with label smoothing.
For DE$\to$EN, MRT increases average OOD BLEU by 0.8 compared to the MLE baseline with label smoothing; for DE$\to$RM the improvement is 0.7 BLEU.
We note that MRT does not consistently improve in-domain performance, which is a first indicator that exposure bias may be more problematic under domain shift.

Our OOD results lag slightly behind those of \citet{mller2019domain}, but note that the techniques employed by them, namely reconstruction \citep{Tu:2017:NMT:3298483.3298684, niu-etal-2019-bi}, subword regularization \citep{kudo-2018-subword}, and noisy channel modelling \citep{li2016mutual} are orthogonal to MRT. We leave the combination of these approaches to future work.

\section{Analysis}

BLEU results indicate that MRT can improve domain robustness.
In this section, we report on additional experiments to establish more direct links between exposure bias and domain robustness, hallucination, and the beam search problem.
Experiments are performed on DE$\to$EN OPUS data.

\subsection{Hallucination}

We manually evaluate the proportion of hallucinated translations on out-of-domain and in-domain test sets.
We follow the definition and evaluation by ~\citet{mller2019domain}, considering a translation a hallucination if it is \textbf{(partially) fluent}, but unrelated in content to the source text (\textbf{inadequate}). We report the proportion of such hallucinations for each system.

\begin{table}[!ht]
\centering
\small
\begin{tabular}{lcc}
\toprule
& \multicolumn{2}{c}{\% hallucinations (BLEU)} \\
\cmidrule{2-3}
system   & out-of-domain & in-domain \\
\midrule
MLE w/o LS  & 35\% \phantom{0}(9.7) & 2\% (58.3) \\
+MRT & 29\% (10.2) & - \\
MLE w/ LS & 33\% (11.2) & 1\% (58.9)\\
+MRT & 26\% (12.0) & - \\
\bottomrule
\end{tabular}
\caption{Proportion of hallucinations and BLEU on out-of-domain and in-domain test sets. DE$\to$EN OPUS.}
\label{tab:3}
\end{table}

Results in Table~\ref{tab:3} confirm that hallucinations are much more pronounced in out-of-domain test sets (33--35\%) than in in-domain test sets (1--2\%).
MRT reduces the proportion of hallucinations on out-of-domain test sets (N=500 for each system; reductions statistically significant at $p<0.05$) and improves BLEU.
Note that the two metrics do not correlate perfectly: MLE with label smoothing has higher BLEU (+1) than MRT based on MLE without label smoothing, but a similar proportion of hallucinations.
This indicates that label smoothing increases translation quality in other aspects, while MRT has a clear effect on the number of hallucinations, reducing it by up to 21\% (relative).

\begin{table}
\begin{center}
\small
\begin{tabular}{ll}
\toprule
source & \textbf{Wir haben} ihn \textbf{gefunden}. \\
reference & \textbf{We found} him.\\
MLE & Do not pass it. \\
MRT & \textbf{We have found} it.  \\
\midrule
source & So höre \textbf{nicht} auf die \textbf{Ableugner}.\\
reference & So hearken \textbf{not} to those who \textbf{deny} (the Truth).\\
MLE & \textbf{Do not} drive or use machines. \\
MRT & \textbf{Do not} apply to \textbf{dleugner}.  \\
\bottomrule
\end{tabular}
\caption{Out-of-domain translation examples. MLE hallucinates in both examples; MRT was rated more adequate in top example, less fluent in bottom one.}
\label{tab:6}
\end{center}
\end{table}

A closer inspection of segments where the MLE system was found to hallucinate shows that some segments were scored higher in adequacy with MRT, others lower in fluency.
One example for each case is shown in Table~\ref{tab:6}.
Even the example where MRT was considered disfluent and inadequate actually shows an attempt to cover the source sentence: the source word `Ableugner' (denier) is mistranslated into `dleugner'.
We consider this preferable to producing a complete hallucination.

\subsection{Uncertainty Analysis}
\label{sec:4.2}
\begin{figure*}[ht]
    \centering
    \includegraphics[width=15cm]{5.pdf}
    \caption{Per-token probability of out-of-domain reference translations and in-domain distractors (first two graphs share legend). Rightmost plot shows direct comparison for MLE baseline and final MRT model. DE$\to$EN OPUS .}
    \label{fig:1}
\end{figure*}

Inspired by~\citet{ott2018analyzing}, we analyse the model's uncertainty by computing the average probability at each time step across a set of sentences.
Besides the reference translations, we also consider a set of `distractor' translations, which are random sentences from the in-domain test set which match the corresponding reference translation in length.

In Figure~\ref{fig:1}, we show out-of-domain results for an MLE model and multiple checkpoints of MRT fine-tuning.
The left two graphs show probabilities for references and distractors, respectively.
The right-most graph shows a direct comparison of probabilities for references and distractors for the MLE baseline and the final MRT model.
The MLE baseline assigns similar probabilities to tokens in the references and the distractors.
Only for the first time steps is there a clear preference for the references over the (mostly random!) distractors.
This shows that error propagation is a big risk: should the model make a wrong prediction initially, this is unlikely to be penalised in later time steps.

MRT tends to increase the model's certainty at later time steps\footnote{The uncertainty of the baseline is due to label smoothing.}, but importantly, the increase is sharper for the reference translations than for the distractors.
The direct comparison shows a widening gap in certainty between the reference and distractor sentences.\footnote{For intermediate checkpoints, see Appendix, Figure~\ref{fig:3}.}
In other words, producing a hallucination will incur a small penalty at each time step (compared to producing the reference), presumably due to a higher reliance on the source signal, lessening the risk of error propagation and hallucinations.

Our analysis shows similar trends on in-domain references.
However, much higher probabilities are assigned to the first few tokens of the references than to the distractors.
Hence, it is much less likely that a hallucination is kept in the beam, or will overtake a good translation in overall probability, reducing the practical impact of the model's over-reliance on its history.\footnote{Figures are shown in the Appendix (Figure~\ref{fig:4}).}

\subsection{Beam Size Analysis}
\label{sec:4.3}

Figure~\ref{fig:1} shows that with MLE, distractor sentences are assigned lower probabilities than the references at the first few time steps, but are assigned similar, potentially even higher probabilities at later time steps.
This establishes a connection between exposure bias and the beam search problem, i.e.\ the problem that increasing the search space can lead to worse model performance.\footnote{The beam search problem has previously been linked to length bias \citep{yang-etal-2018-breaking,murray-chiang-2018-correcting} and the copy mode ~\citep{ott2018analyzing}. We consider hallucinations another result of using large search spaces with MLE models.} With larger beam size, it is more likely that hallucinations survive pruning at the first few time steps, and with high probabilities assigned to them at later time steps, there is a chance that they become the top-scoring translation.

We investigate whether the beam search problem is mitigated by MRT.
In Table~\ref{tab:beam}, we report OOD BLEU and the proportion of hallucinations with beam sizes of 1, 4 and 50.
While MRT does not eliminate the beam search problem, performance drops less steeply as beam size increases.
With beam size 4, our MRT models outperform the MLE baseline by 0.5-0.8 BLEU; with beam size 50, this difference grows to 0.6-1.5 BLEU.
Our manual evaluation (N=200 for each system for beam size 1 and 50) shows that the proportion of hallucinations increases with beam size, and that MRT consistently reduces the proportion by 11-21\% (relative).
For the system with label smoothing, the relative increase in hallucinations with increasing beam size is also smaller with MRT (+33\%) than with MLE (+44\%).

\begin{table}[ht]
\centering
\small
\begin{tabular}{lccc}
\toprule
& \multicolumn{3}{c}{BLEU (\% hallucinations)} \\
\cmidrule{2-4}
system   & $k=1$ & $k=4$ & $k=50$ \\
\midrule
MLE w/o LS & \phantom{0}8.9 (28\%) & \phantom{0}9.7 (35\%) & \phantom{0}9.3 (37\%)\\
+MRT & \phantom{0}9.1 (24\%) & 10.2 (29\%) & \phantom{0}9.9 (33\%)\\
MLE w/ LS & 10.6 (27\%) & 11.2 (33\%) & \phantom{0}9.4 (39\%)\\
+MRT & 11.3 (24\%) & 12.0 (26\%) & 10.9 (32\%)\\
\bottomrule
\end{tabular}
\caption{Average OOD BLEU and proportion of hallucinations with different beam sizes $k$. DE$\to$EN OPUS.}
\label{tab:beam}
\end{table}

\section{Conclusions}

Our results and analysis show a connection between the exposure bias due to MLE training with teacher forcing and several well-known problems in neural machine translation, namely poor performance under domain shift, hallucinated translations, and deteriorating performance with increasing beam size.
We find that Minimum Risk Training, which does not suffer from exposure bias, can be useful even when it does not increase performance on an in-domain test set: it increases performance under domain shift, reduces the number of hallucinations substantially, and makes beam search with large beams more stable.

Our findings are pertinent to the academic debate how big of a problem exposure bias is in practice -- we find that this can vary substantially depending on the dataset --, and they provide a new justification for sequence-level training objectives that reduce or eliminate exposure bias.
Furthermore, we believe that a better understanding of the links between exposure bias and well-known translation problems will help practitioners decide when sequence-level training objectives are especially promising, for example in settings where the test domain is unknown, or where hallucinations are a common problem.

\section*{Acknowledgments}

Chaojun Wang was supported by the UK Engineering and Physical Sciences Research Council (EPSRC) fellowship grant EP/S001271/1 (MTStretch).
Rico Sennrich acknowledges support of the Swiss National Science Foundation (MUTAMUR; no.\ 176727).
This project has received support from Samsung Electronics Polska sp. z o.o. - Samsung R\&D Institute Poland.

\clearpage


text:
\section{Introduction and related works \label{SectionIntroduction}}

Reliable computation with unreliable devices, or in the presence of noise, has been the subject of numerous works within the vast field of fault-tolerant computing. Computation can be made reliable using information and  component/gate redundancy. Some works aim at identifying theoretical "boundaries" on the amount of necessary and/or sufficient redundancy to achieve reliability.   Recent references (see for example, \cite{Spielman96}, \cite{Gacs05},\cite{Hadjicostisverghese05}, \cite{RachlinSavage08}) continue to extend the stream opened by Von Neumann's seminal paper \cite{VonNeumann56}. These works identify bounds (e.g., depth and size of circuits) and  propose frameworks to design reliable computations mainly thanks  to gate redundancy.  These papers do not address the question of boundaries about  information redundancy. This subject, through the concept of capacity and coding theorems, has been thoroughly studied for data communication. It has not been the case for computation although the two problems (noisy computation and noisy communication) are very close: the matter is to retrieve expected values from the outputs of a random process.

The question of whether a noisy computation possesses a capacity (or equivalently whether some coding theorems for noisy computation hold and in which cases) has been raised by P. Elias in 1958  \cite{Elias58}. There is a practical consequence in answering positively this question. This would  mean that, given a noisy implementation of an expected (e.g., Turing computable) function, it is possible to find families of efficient  input codes which asymptotically allow an almost perfect computation. Efficiency means having an input encoding rate which could remain strictly positive or arbitrarily close to a capacity when the length of the code tends to infinity.

In a strongly constrained context (independent encoding of operands for bit-by-bit boolean operations), P. Elias obtained negative first results on the existence of a noisy computation capacity (\cite{Elias58}) . This work was deepened by Peterson and Rabin in \cite{PetersonRabin59} and by Winograd in \cite{Winograd62}. One of the major conclusions of these studies was that reliable computation with positive rate (the ratio $\frac{k}{n}$ of encoding $k$-length input block in $n$-length blocks of binary symbols) in the presence of noise  is not possible for some boolean operations (e.g., AND) under some assumptions (independent coding of operands, bijective decoding and bit-by-bit operation). This led to the conclusion that, under these assumptions, there is no capacity for such noisy operations. It is worth noting that the assumptions were made to forbid the reliable encoder and decoder to "participate" to the computation of the expected operation.

In \cite{Ahlswede84}, Ahlswede went into the subject in greater depth and made an important contribution  in characterizing  contexts  in which  a capacity for noisy computations cannot exist. It appears that the characteristics of the decoding function play a fundamental role. If the inverse of the decoding function is injective and monotonic  then weak converse theorems hold for the average and maximal error probabilities. If, in addition, the inverse of the decoding function preserves the logical AND (this implies monotonicity), then strong converses hold. These theorems state that the rate of encoding tends to $0$ when the block code length tends to infinity. The hypotheses made in \cite{Elias58,PetersonRabin59,Winograd62}, i.e., independent encoding of operands and bijective decoding, imply monotonicity of the inverse of the decoding function. On this aspect, \cite{Ahlswede84} supersedes \cite{Elias58,PetersonRabin59,Winograd62}.

Nevertheless, these negative results do not imply the absolute impossibility to identify a capacity for noisy computation. They characterize   codes, encoding and decoding processes which cannot open this ability. To  define a capacity for noisy computations, assumptions must be relaxed.

To the author's best knowledge, the first positive answer given through a definition of a capacity of a noisy computation (in fact similar to the one for a noisy channel) and a coding theorem came from  Winograd and Cowan in \cite{WinogradCowan63}.   In \cite{WinogradCowan63}, the entropy $H(X|F(X))$ of the input source conditioned by the noisy computation output is assessed as a noise measure. As it is the equivocation between the {\em noisy output} and the {\em input}, this quantity is not relevant, in full generality,  as the equivocation due to the sole noise:  it encompasses also the amount of information lost by computation. But, in a special case of noisy functions called decomposable modules, $H(X|F(X))$ actually measures the equivocation due to noise.  Decomposable modules are noisy functions  which  can be modeled by a perfect function followed by a noisy communication channel: the error probability depends on the desired output value rather than on the input value. These peculiar  functions, though noisy, make the context equivalent to  that where the reliable encoder computes and encodes the expected function result before communication through a noisy channel. Due to the restriction of considering decomposable modules, \cite{WinogradCowan63} did not completely succeed in proposing a noisy computation capacity  in a general scope (\cite{WinogradCowan63}, theorem 6.3, pages 47-48 ).

Noisy computation capacity is also considered in reliable reconstruction of a function of sources over a multiple access channel. Much more recently, a definition of noisy computation capacity is established by Nazer and Gastpar in \cite{NazerGastpar07} and is totally consistent with the one proposed here.  Nazer and Gastpar demonstrate the possible advantages of joint source-channel coding of multiple sources over a separation-based scheme, allowing a decoder to retrieve a value which is a function of input sources. This context makes relevant the proposed distributed encoding process which perfectly performs a computation equivalent to the desired function. The encoder outputs are then transmitted through a noisy MAC to a decoder (see proofs of Theorems 1 and 2 of \cite{NazerGastpar07}). This also models a noisy computation as a perfect computation followed by a noisy transmission of the result. It can be noticed that \cite{WinogradCowan63} and \cite{NazerGastpar07} relax the assumptions of \cite{Elias58, PetersonRabin59, Winograd62, Ahlswede84} in a similar way: all goes as if the operands are jointly coded into an encoded form of the expected function result before being handled by a noisy communication channel.

The present  paper establishes a model setting down the problem of noisy computation (section \ref{SectionNoisyComputation}), a definition of the capacity of a noisy computation with respect to an expected function and a coding lemma (section \ref{SectionCapacityofaNoisyComputation}).  A  model for reliable computation is  given, section \ref{SectionReliableComputation}. Based on this model, a joint source-computation coding theorem and its converse are stated and proved in Section \ref{SectionCodingTheorem}. This theorem aims at formally capturing  practical approaches in which reliable computation of a function $g$ is obtained thanks to a noisy apparatus $F$ computing with noise a function $f$ (e.g., a regular arithmetic addition $g$ obtained from the noisy actual circuit $F$ implementing $f$ which is an addition acting on residue encoded operands, \cite{RaoFujiwara89}). The input source, $f$, $g$, $F$ and the cascade $f^{-1}F$ are supposed to be  AMS  and ergodic one-sided random processes or channels, extending \cite{Simon10} to more general random processes and algorithms. The perfect function $f$ is assumed  {\em unary} (as is a Turing computable function). $n$-ary functions can be modeled as unary ones by concatenating $n$ input values in one "meta"-input and thus modeling a {\em joint coding} of operands. This relaxes the assumptions of  \cite{Elias58, PetersonRabin59, Winograd62, Ahlswede84}.

\section{Model for Noisy Computation \label{SectionNoisyComputation}}

 In this section, the notations used follow \cite{Gray11} and cover countable alphabets, assumed standard (thus conditional probabilities are regular).

Let $X\equiv \{ X_i ; i\in \mathcal{I} \}$ a random process with values in  $(A^\mathcal{I}, \mathcal{B}_{A^\mathcal{I}})$ where $A$ is a countable alphabet and $\mathcal{I}$ a countable set of indexes (e.g., $\mathbb{N}$). $\mathcal{B}_{A^\mathcal{I}}$ denotes the $\sigma$-field generated by the rectangles (chapter 1, \cite{Gray11}). The random process $X$ is the source of inputs of the noisy computation.

The noisy computing device is modeled as a random channel, i.e. a set of conditional probabilities $F \equiv \{ F_x, x \in A^\mathcal{I}\}$, taking $X$ as input and producing as an output a random process $Z\equiv \{ Z_i ; i\in \mathcal{I} \}$  on $(C^\mathcal{I}, \mathcal{B}_{C^\mathcal{I}})$ where $C$ is a countable alphabet. The hookup $P_{XZ}\equiv P_XF$ is the probability measure characterizing the actual noisy computation with an input flow represented by $X$.  From \cite{Gray11}, chapter 2, $\forall O \in \mathcal{B}_{A^\mathcal{I} \times C^\mathcal{I}}$:
$$
P_{XZ} (O) =  \int_{A^\mathcal{I}} P_{Z|X}(O_x|x) dP_X
	 = \int_{C^\mathcal{I}} P_{X|Z}(O_z|z) dP_Z
$$
where $O_x=\left\{z\in C^\mathcal{I} / (x,z)\in O \right\}$. The probabilities $\left\{   P_{Z|X}(. |x), x\in A^\mathcal{I}\right\}$ defines the channel $X \rightarrow Z$ ($F_x \equiv P_{Z|X}(. |x)$) and $\left\{   P_{X|Z}(. |z), z\in C^\mathcal{I}\right\}$ the  "reverse" channel $F^{-1}$.

$P_{XZ}$ fully determines the channels $F$ and $F^{-1}$. Conversely, if $X$ and a set of conditional probabilities $\left\{   P_{Z|X}(. |x), x\in A^\mathcal{I}\right\}$ (i.e. $F$) are given, then $P_{XZ}$ and the output process $Z$ are well defined. A functional notation $Z=F(X)$ will be used below.

The desired (i.e. perfect) computation will be represented by a {\em measurable} function $f: A^\mathcal{I} \rightarrow B^\mathcal{I}$ where $B$ is a countable alphabet.  $Y=f(X)$ is a random process of  distribution $P_Y=P_Xf^{-1}$. The function $f$ defines a deterministic channel $X \rightarrow Y=f(X)$ which is a set of conditional probabilities $\left\{   P_{f(X)|X}(. |x), x\in A^\mathcal{I}\right\}$ (See \cite{Gray11}, chap. 2):
\begin{equation}
	\forall G \in \mathcal{B}_{B^\mathcal{I}},  P_{f(X)|X}(G|x) = 1_{f^{-1}(G)}(x) \text{ 	}P_X\text { a.e.}  \label{EquationDeterministicChannel}
\end{equation}

$F$ and $f$ determine a channel $f(X)\rightarrow Z$ which is a cascade of the reverse channel $f^{-1}$ followed by $F$ (figure~(\ref{FigNoisyComputationModel}).
\begin{figure}[htbp]
	\centering
	\xymatrix{
	\left( A^{\mathcal{I}} , \mathcal{B}_{A^\mathcal{I}}, P_X \right) 	\ar @/^1.5pc/ [rr]^f \ar [dr]_F	&						& \left( B^{\mathcal{I}} , \mathcal{B}_{B^\mathcal{I}}, P_Y  \right) \ar @/_-1pc/  [ll]^{f^{-1}} \ar [dl]^{f^{-1}F}\\
 													&  \left( C^{\mathcal{I}} , \mathcal{B}_{C^\mathcal{I}}, P_Z \right)
 	 }
	\caption{Model for Noisy Computation}
	 \label{FigNoisyComputationModel}
\end{figure}
The noisy computation model should not be understood as a cascade made of a perfect function followed by a noisy channel (as done in \cite{WinogradCowan63} for instance). The cascade $f^{-1}F$ is an "artifact" on which  a channel coding theorem will be invoked to build a code for the noisy computation.

From the expression giving the probabilities of a cascade of channels  (\cite{Gray11}, chap. 2) and applying (\ref{EquationDeterministicChannel}), it easily comes that $P_{Y|Z}\equiv P_{X|Z}f^{-1}$.

 If $f$ is bijective and $A=B$, the model is the noisy channel one.

\section{Capacity of a Noisy Computation \label{SectionCapacityofaNoisyComputation}}

The alphabets are now assumed finite and the $\sigma$-fields $\mathcal{B}_A$, $\mathcal{B}_B$  and $\mathcal{B}_C$ are the sets of subsets of $A$, $B$ and $C$. $f^n(X^n)$  stands for the random value $Y^n$ (the $n$ first symbols of $Y$), $F^n(X^n)$ for $Z^n$. $[A\times B, Xf ; A\times C, XF]$ denotes the noisy computation $F$ of $f$ on $X$.
{\em $X$, $F$, $f$ and the cascade $f^{-1}F$ are assumed AMS and ergodic}. Thus entropy rates are limits and ergodic theorems hold (\cite{Gray11}).

The Feinstein's theorem is reminded (\cite{Gray11}, chap. 12):
\begin{theorem}[Feinstein's theorem]
Let $[A\times B, \mu\nu]$ be an AMS and ergodic hookup of a source $\mu$ and channel $\nu$. Let $\overline{I}_{\mu\nu}=\overline{I}_{\mu\nu}(X;Y)$ denote the average mutual information rate and assume that $\overline{I}_{\mu\nu}$ is finite (this is the case if the alphabets are finite.). Then for any $R<\overline{I}_{\mu\nu}$ and for any $\epsilon >0$, there exists, for $n$ large enough, a code $\left\{ (\omega_i, \Gamma_i) \in A^n \times \mathcal{B}_{B^n} , i=1,\cdots,M \right\}$
such that $M = \lfloor e^{nR} \rfloor $ and $\forall i=1,\cdots,M,  \hat{\nu}^n(\Gamma_i^c|\omega_i) \leq \epsilon$
\end{theorem}
$\hat{\nu}^n$ is the channel induced by the source $\mu$, (\cite{Gray11}, chap. 12).

The following definition introduces the {\em typical input rate} which will be shown to be the rate at which a source should produce typical inputs for a random process to allow to recover the desired function results by decoding.

\begin{definition}
The {\em typical input rate} of the source $X$ for the noisy computation $F$ with respect to the perfect function $f$ is the following limit, denoted $\overline{B}(X,f,F)$:
$$
	lim_{n\rightarrow \infty }\frac{ H(X^n|f^n(X^n)) + I(F^n(X^n);f^n(X^n))}{n} \nonumber
$$
\end{definition}
Taking into account the hypothesis ( AMS and ergodic source and channels), $\overline{B}(X,f,F)$  is  well defined. Simple algebra gives $\overline{B}(X,f,F)=\overline{H}(X) - \overline{H}(f(X)|F(X))$

\begin{definition} \label{DefinitionFeinsteinCodeForNoisyComputation}
Let $[A\times B, Xf ; A\times C , XF]$ be a noisy computation on finite alphabets $A$, $B$ and $C$. A {\em $[M,n,\epsilon]$-Feinstein code} for the noisy computation $[A\times B, Xf ; A\times C , XF]$ is a set $\left\{(A^n_i, \Gamma^n_i)\in \mathcal{B}_{A^n}\times \mathcal{B}_{C^n}, i=1,\cdots,M \right\}$ such that:
\begin{enumerate}
	\item $\hat{P}^n_{F(X)|X}(\Gamma_i^c|x^n) \leq \epsilon$ for any $x^n\in A^n_i$, $i=1,\cdots,M$
	\item $\forall i=1,\cdots,M, \exists y_i \in B^n \text{ such that  } A^n_i =  (f^n)^{-1}(y_i)$
\end{enumerate}
$F^n(X^n)$ {\em $\epsilon$-reliably computes} $f^n(X^n)$ on the code $\left\{(A^n_i, \Gamma^n_i), i=1,\cdots,M \right\}$.
\end{definition}

\begin{lemma} [Feinstein's theorem for Noisy Computation] \label{TheoremNoisyComputation}
Let $[A\times B, Xf ; A\times C , XF]$ be a noisy computation on finite alphabets. For any $R<\overline{B}(X,f,F)$, for any $\epsilon>0$, n large enough, there exists a $[\lfloor e^{n(R-\overline{H}(X|f(X)))} \rfloor,n,\epsilon]$-Feinstein code for $F^n(X^n)$ to $\epsilon$-reliably computes $f^n(X^n)$.
\end{lemma}
\begin{IEEEproof}
$R<\overline{B}(X,f,F) \Rightarrow R'=R-\overline{H}(X|f(X)) < \overline{I}(f(X),F(X))$. Then, thanks to the Feinstein's theorem, since $P_{f(X),F(X)}$ is AMS and ergodic (by assumption), for $n$ large enough, there exists a Feinstein code $\{ (y_i, \Gamma_i) \in B^n \times \mathcal{B}_{C^n} ; i=1,\cdots, M\}$ such that $M=\lfloor e^{nR'} \rfloor$ and $\forall i=1,\cdots,M ; \hat{P}^n_{F(X)|X}(\Gamma_i^c|y_i) \leq \epsilon$

Let $x$ belong to $ (f^n)^{-1}(y)$. Considering  the cascade $X^n\rightarrow f^n(X^n)\rightarrow F^n(X^n)$, for any $k=1,\cdots,M$:
$$
\hat{P}^n_{F(X)|X}(\Gamma_k^{c}|x) = \int_{B^{\mathcal{I}}} \hat{P}^n_{F(X)|f(X)}(\Gamma_k^{c}|\underline{y}) d\hat{P}^n_{f(X)|X}(\underline{y}|x)
$$
\begin{multline}
 \hat{P}^n_{F(X)|X}(\Gamma_k^{c}|x)  =
 	\int_{\left\{y_k \right \}} \hat{P}^n_{F(X)|f(X)}(\Gamma_k^{c}|\underline{y}) d\hat{P}^n_{f(X)|X}(\underline{y}|x) \nonumber \\
				    + \int_{\left\{y_k \right \}^c} \hat{P}^n_{F(X)|f(X)}(\Gamma_k^{c}|\underline{y}) d\hat{P}_{f(X)|X}(\underline{y}|x) \nonumber \\
			         		 \leq   \hat{P}^n_{f(X)|X}\left(  \left\{y_k \right \} |x \right).\epsilon +\hat{P}^n_{f(X)|X}\left(   \left\{y_k \right \}^c |x \right) \nonumber
\end{multline}
If $x\in (f^n)^{-1}(y_k)$ then $\hat{P}^n_{f(X)|X}\left(  \left\{y_k \right \} |x \right)=1$ and $\hat{P}^n_{f(X)|X}\left(  \left\{y_k \right \}^c |x \right)=0$, hence $\forall x\in f^{-1}(y_k)$,  $ \hat{P}^n_{F(X)|X}\left( \Gamma_k^c | x_k  \right) \leq \epsilon$
\end{IEEEproof}

We can conclude this section by the definition of the {\em typical input capacity} of a noisy computation.

\begin{definition}
The typical input capacity of the noisy function $F$ with respect to the perfect function $f$ is $C_f(F) = \sup_{\text{AMS erg  } P_X} \overline{B}(X,f,F)$, the supremum is over all AMS and ergodic sources $X$.
\end{definition}

The equivalent expression $C_f(F) = \sup_{\text{AMS erg  } P_X} \left[ \overline{H}(X) - \overline{H}(f(X)|F(X))\right]$  shows that this capacity boils down to the "usual" channel capacity when $f$ is a bijection, in which case $\overline{H}(f(X)|F(X))=\overline{H}(X|F(X))$.

\section{Reliable computation\label{SectionReliableComputation}}
There is a need (\cite{Elias58,Winograd62}) to constrain the encoding and decoding processes to avoid  the following cases:
\begin{itemize}
	\item either an (assumed perfect) encoder which computes the expected function, encodes the result before transmission through the random process (considered as a noisy transmission channel)
	\item or an encoder which encodes input values for reliable transmission through the random process (considered here also as a noisy transmission channel) and a decoder (assumed also reliable) which decodes (almost perfectly) and computes (perfectly) the expected function.
\end{itemize}

Considering that a computation $g$  is a "true" computation if the entropy is reduced ($\overline{H}(X'|g(X'))>0 \Rightarrow \overline{H}(g(X')) < \overline{H}(X')$ (else it is communication), the model must be targeted to be  mainly relevant for non-injective functions (i.e. $\overline{H}(X'|g(X'))>0$). For injective $g$, this becomes  the classical reliable transmission model. With the constraint that both encoding and decoding are based on injections (in a sense made precise below) then the encoder and the decoder cannot compute (at least totally) the desired function as they do not reduce entropy.

The proposed model of the complete process to reliably compute a  function $g :A'^{\mathcal{I}} \rightarrow B'^{\mathcal{I}}$ acting on a source $X'$, thanks to a  noisy implementation $F$ of a  function $f: A^{\mathcal{I}} \rightarrow C^{\mathcal{I}}$ is the following:
\begin{itemize}
	\item {\bf encoding}: let $X^n$ be the $n^{th}$ extension of an  source for which we have a maximal code $(A^n_i,\Gamma_i)_{i=1,\cdots,M}$ allowing to $\epsilon$-reliably compute $f^n(X^n)$ by $F^n(X^n)$ (cf lemma~\ref{TheoremNoisyComputation} and definition~\ref{DefinitionFeinsteinCodeForNoisyComputation}) ; a typical $k$-sequence $x'$ of $X'^k$ is encoded into a {\em typical given $y_i$} (this important assumption is discussed in the conclusion) $n$-sequence of $X^n$ by a injective function, say $\mathcal{U}$, such that $\mathcal{U}(x') \in A_i^n$ for some $i=1,\cdots,M$
	\item {\bf computation of the noisy function}: $F^n$ is applied to $\mathcal{U}(x')$ producing a typical $n$-sequence $F^n(\mathcal{U}(x'))$ of $F^n(X^n)$ where $F^n(\mathcal{U}(x'))$ belongs to a given $\Gamma_i$ (with probability greater than $1-\epsilon$)
	\item {\bf decoding}: the first step is to associate to $F^n(\mathcal{U}(x'))$ the typical $n$-sequence $y_i$ of $f^n(X^n)$ corresponding to $\Gamma_i$, the second step is to apply to $y_i$ a function $\mathcal{V} : \{\mathbf{y_1},\ldots,\mathbf{y_M}\} \rightarrow \{\text{typical k-sequences of }g^k(X'^k)\} $ such that $\mathcal{V}(y_i)=g^k(x')$
\end{itemize}
A decoding error occurs when one obtains a $n$-sequence $y_j$ (or equivalently a  $\Gamma_j$) such that $\widehat{g^k(x')} =  \mathcal{V}(y_j) \neq g^k(x')$

To be able to define a decoding function $\mathcal{V}$ (i.e, a deterministic decoding), the encoding function $\mathcal{U}$ has to be such that the typical (given $y_i$) $n$-sequences of one $A_i^n = (f^n)^{-1} (y_i)$ ($ y_i \in \{y_1,\ldots,y_M\}$) are used for encoding typical $k$-sequences of {\em only one} $(g^k)^{-1} (z)$, $z$ typical $k$-sequence of $g^k(X'^k)$.

We also require that $\mathcal{V}$ be an injection (as we have required from $\mathcal{U}$).

The typical $k$-sequences of a $(g^k)^{-1} (z)$, $z$ typical $k$-sequence of $g^k(X'^k)$, are encoded in typical (given $y_i$) $n$-sequences of one and only one $A_i^n$ . So, if $x'_1$ and $x'_2$ are two typical $k$-sequences of $X'^k$:
$$
f^n(\mathcal{U}(x'_1)) = f^n(\mathcal{U}(x'_2)) \Leftrightarrow g^k(x'_1) = g^k(x'_2)
$$

The model fulfills the constraints identified above. The encoder implements an injection and thus cannot compute the desired function $f$ nor $g$ (if $f$ and $g$ are not injective). The same comment applies to the injective decoding step  $\mathcal{V}$.

\section{A coding theorem and its converse\label{SectionCodingTheorem}}

The sources, functions, noisy function and the cascade $f^{-1}F$ are assumed AMS and ergodic.

\begin{definition}   \label{DefinitionEncodingRate}
With the notations of section \ref{SectionReliableComputation}, the ratio $R=\frac{k.H(X')}{n}$ is called the {\em typical encoding input rate}. A rate $R$ is said to be {\em achievable with respect to the function $f$} if there exists a sequence of codes of size $n$ such that the maximal probability of decoding error tends to 0 as $n$ tends to infinity.
\end{definition}

\begin{theorem} \label{TheoremCoding}
If $R < C_f(F)$, then $R$ is achievable w.r.t $f$.
\end{theorem}

\begin{IEEEproof}

This proof, although identical to that in \cite{Simon10}, is given as it includes the starting point for the proof of the converse theorem.
 First, it is shown that the injective encoding of typical $k$-sequences of a set $(g^k)^{-1} (z)$ on typical (given $y_i$) $n$-sequences belonging to $A_i^n$ is possible for suitably chosen $k$ and $n$ (lossless coding). Secondly, it is shown that, at encoding input rates below capacity and for $k$ and $n$ suitably chosen,  the sets $A_i^n$ are almost as many as the sets $(g^k)^{-1} (z)$.

Let $\delta'' >0$. Since $\mathbb{Q}$ is dense in $\mathbb{R}$, there exist $k$ and $n$ such that:
$$
\frac{\overline{H}(X'|g(X'))}{\overline{H}(X|f(X))} < \frac{n}{k} <\frac{\overline{H}(X'|g(X'))+\delta''}{\overline{H}(X|f(X))}
$$
Moreover, $k$ and $n$ can be chosen as large as needed. Thus:
\begin{equation}
\frac{k.\overline{H}(X'|g('X))}{H(X|f(X))} < n <\frac{k.(\overline{H}(X'|g(X'))+\delta'')}{\overline{H}(X|f(X))}  \label{InequalitySourceCoding}
\end{equation}
We can choose $\delta, \delta' >0$ and $0<\epsilon<1/2$ small enough for:
\begin{multline*}
	\frac{k.\left( \overline{H}(X'|g(X'))+\delta \right) - log(1-2\epsilon)}{\overline{H}(X|f(X))-\delta'} < n \\
	<\frac{k.(\overline{H}(X'|g(X'))+\delta+\delta'')}{\overline{H}(X|f(X))+\delta'}
\end{multline*}
giving
\begin{multline*}
	k.(\overline{H}(X'|g(X'))+\delta)  < log(1-2\epsilon) + n.(\overline{H}(X|f(X))-\delta') \\
 	< n.(H(X|f(X))+\delta') < k.(H(X'|g(X'))+\delta+\delta'')
\end{multline*}
If $\nu_1$ is the number of typical $k$-sequences of $(g^k)^{-1} (z)$ and $\nu_2$ is the number of typical (given $y_i$) $n$-sequences in an $A_i^n$, we have (by conditional AEP):
\begin{multline*}
	\nu_1 < e^{k.(\overline{H}(X'|g(X'))+\delta)}  <(1-2\epsilon)e^{ n.(\overline{H}(X|f(X))-\delta')} \\
	<\nu_2 <e^{ n.(\overline{H}(X|f(X))+\delta')} <e^{ k.(\overline{H}(X'|g(X'))+\delta+\delta'')}
\end{multline*}
It is thus possible to find an injection from the set of typical $k$-sequences of $(g^k)^{-1} (z)$ on the subset of typical sequences (given $y_i$) of $A_i^n$. This shows the first step.

Assume that $R=k\overline{H}(X')/n< \overline{H}(X)-\overline{H}(f(X)|F(X)) \leq C_f(F)$. Such a  $X$ exists by definition of $C_f(F)$.  So
\begin{multline*}
	k(\overline{H}(g(X'))+\overline{H}(X'|g(X'))) < \\
 	 n.(\overline{H}(f(X))-\overline{H}(f(X)|F(X)))+n.\overline{H}(X|f(X))
\end{multline*}

By (\ref{InequalitySourceCoding}), $n.\overline{H}(X|f(X)) - k.\overline{H}(X'|g(X')) < k.\delta''$
thus
$$
k\overline{H}(g(X')) < n.(\overline{H}(f(X))-\overline{H}(f(X)|F(X)))+ k.\delta''
$$
$\epsilon_1,\delta''' >0$ can be chosen small enough in order to get:
$$
e^{k(\overline{H}(g(X'))+\delta''')} < e^{n.(\overline{H}(f(X))-\overline{H}(f(X)|F(X))+ \frac{k}{n}.\delta''-\epsilon_1)}
$$

If $\nu_3$ is the number of typical $k$-sequences of $g^k (X'^k)$ and $M$ is the size of the code (i.e., the number of $(A_i^n,\Gamma_i)$), we have (by AEP and Lemma~\ref{TheoremNoisyComputation}):
\begin{multline*}
\nu_3 < e^{k(\overline{H}(g(X'))+\delta''')} \\
	< e^{n.(\overline{H}(f(X))-\overline{H}(f(X)|F(X))+ \frac{k}{n}.\delta''-\epsilon_1)} <M
\end{multline*}
\end{IEEEproof}

The assumed model, by the constraints on encoding, implies that the best ratio (i.e., the smaller) $\frac{n}{k}$ of encoding respects the inequality~(\ref {InequalitySourceCoding}):  $\frac{\overline{H}(X'|g(X'))}{\overline{H}(X|f(X))} < \frac{n}{k}<\frac{\overline{H}(X'|g(X'))+\delta''}{\overline{H}(X|f(X))}$. Let $\gamma=\frac{\overline{H}(X'|g(X'))}{\overline{H}(X|f(X))}$. To respect the encoding constraints (typical sequences are "injectively" encoded into typical sequences), a rate $R=\frac{k}{n}\overline{H}(X')$ must be such that $R\leq \frac{\overline{H}(X')}{\gamma}$

\begin{theorem} \label{TheoremConverse}
If $R>C_f(F)$, there is no code such that the error probability tends to $0$ as $n \rightarrow \infty$
\end{theorem}

\begin{IEEEproof}
The decoding is deterministic then:
\begin{enumerate}
	\item $F^n(X^n)\rightarrow f^n(X^n)  \rightarrow g^k(X'^k)$ is a Markov Chain  thus $g^k(X'^k)\rightarrow f^n(X^n)  \rightarrow F^n(X^n)$ is a Markov Chain
	\item $f^n(X^n)\rightarrow F^n(X^n)  \rightarrow \widehat{g^k(X'^k)}$ is a Markov Chain
\end{enumerate}
Hence $g^k(X'^k)\rightarrow f^n(X^n)  \rightarrow F^n(X^n)  \rightarrow \widehat{g^k(X'^k)} $ is a Markov Chain. This implies that, $\forall n,k \text{  such that }\frac{k}{n} \leq \frac{1}{\gamma}$:
\begin{multline}
I(g^k(X'^k) ; \widehat{g^k(X'^k)})  \leq  I(f^n(X^n) ; F^n(X^n)) \text{ hence}\nonumber \\
 H(g^k(X'^k)) - H(g^k(X'^k)|\widehat{g^k(X'^k)})  \leq  \nonumber \\
                                           			   I(f^n(X^n) ; F^n(X^n)) \nonumber \\
\Rightarrow H(X'^k)- H(X'^k | g^k(X'^k)) - H(g^k(X'^k)|\widehat{g^k(X'^k)})  \leq \nonumber \\
                I(f^n(X^n) ; F^n(X^n)) \nonumber
\end{multline}
by Fano's inequality:
\begin{multline}
H(X'^k)- H(X'^k | g^k(X'^k)) \nonumber \\
	- (H_2(P_e(k))+k.P_e(k).log(|B'|))  \nonumber \\
		\leq  I(f^n(X^n) ; F^n(X^n)) \nonumber \\
\Rightarrow \frac{H(X'^k)}{k}  - \frac{H_2(P_e(k))}{k} - P_e(k).log(|B'|)  \leq  \nonumber \\
	\frac{H(X'^k | g^k(X'^k))}{k} + \frac{ I(f^n(X^n) ; F^n(X^n))}{k}\nonumber \\
\Rightarrow \frac{H(X'^k)}{k}  - \frac{H_2(P_e(k))}{k} - P_e(k).log(|B'|)  \leq \nonumber \\
	 \frac{H(X'^k | g^k(X'^k))}{k} + \gamma . \frac{ I(f^n(X^n) ; F^n(X^n))}{n}\nonumber
\end{multline}

If the error probability is asymptotically $0$ (i.e., $lim_{k \rightarrow \infty}(P_e(k)=0$) then necessarily (letting $k$ and $n$ tend to infinity): $\overline{H}(X')  \leq  \overline{H}(X' | g(X')) + \gamma . \overline{I}(f(X) ; F(X))$.
But $\overline{H}(X' | g(X')) = \gamma \overline{H}(X | f(X))$, then:
$$
\overline{H}(X')  \leq  \gamma .\left( \overline{H}(X | f(X)) + \overline{I}(f(X) ; F(X)) \right)
$$
since $R\leq \frac{\overline{H}(X')}{\gamma}$, we obtain
$$
R \leq \overline{H}(X | f(X)) + \overline{I}(f(X) ; F(X)) \leq C_f(F)
$$
Thus if $R > C_f(F)$ then the error probability does not vanish.
\end{IEEEproof}

\section{Discussion and conclusion \label{SectionConclusion}}

The coding lemma, the coding theorem and its converse assume that the sources ($X$ and $X'$),  the channels ($F$, $f$, $g$) and the cascade $f^{-1}F$ are AMS and ergodic. Cases can be identified where it is possible to derive the AMS property and ergodicity of  $f^{-1}F$ from properties of $X$, $F$ and $f$ (e.g., if $X$, $F$ are stationary and weakly mixing and $f$ AMS and ergodic then $f^{-1}F$ is AMS and ergodic). Due to lack of space, this question is not addressed here, neither the identification of classes of AMS and ergodic functions $f$ and $g$.

The model of reliable computation assumes that the encoder and the decoder are perfectly reliable. This assumption could be justified by quoting from \cite{Winograd62} {\em "The computation system [model] was devised for the sole purpose of studying the relation of information theory of reliable automata"}. Moreover we could argue in addition that if the complexity of the computation device is of a much greater magnitude than that of the encoder and decoder then the unreliability of the encoder and decoder have almost no impact on the overall reliability of the computation and thus can be neglected. For complex systems (e.g., based on significant software volume), this is quite realistic. In any case, it is impossible to overcome the fact that the reliability reached is at the best the reliability of the final decoding device. The only way is to built a intrinsically reliable enough decoder (for example thanks to gate redundancy).  A noisy encoder is a noisy computation itself and thus can be handled from the point of view of "cascaded noisy computations". This is outside of the scope of the present paper.

The proposed model of reliable computation involves two perfect functions $ g$ and $f$. This is intended to capture major real cases as already mentioned. Another motivation to use an "ancillary" function $f$ in the model is that this is an efficient way to define an input code,  meaning a family of subsets $(A^n_i)_{i=1,\cdots,M}$, that do not overlap and whose "images" by the noisy function do not overlap "too much" (i.e., fall into disjoint $\Gamma_i$ with high probability). Defining such family is defining (partially) a function $f$ by picking, for each $i$ an $y_i$ and stating $f^{-1}(y_i)=A_i$. In addition $f$ allows a characterization of a kind of size of the sets $(f^n)^{-1}(y_i)$ through the conditional entropy rate $\overline{H}(X|f(X))$. This motivates also  the constraint of coding by conditionally typical sequences. While the sets $(f^n)^{-1}(y_i)$ are "balanced", for large $n$, with respect to the number of (conditionally) typical sequences they contain, their cardinalities might be very different and bounds are not straightforward to obtain. Thus, the use of all possible  elements of $(f^n)^{-1}(y_i)$ forbids to characterize all the $(f^n)^{-1}(y_i)$ by the same number measuring the "encoding" power. The same difficulty forbids to state a converse as well. The "encoding by conditionally typical sequences" trick overcomes this difficulty.

\bibliographystyle{IEEEtran}

redpajama_set_name: RedPajamaArXiv


text:
\section{Introduction}

The discovery of FeAs-based~\cite{Kamihara2008} and FeSe-based~\cite{Hsu2008} superconductors created an exciting platform for better understanding the physics of high-temperature superconductivity. The observation of wide-range of transition temperature ($T_\text{c}$), with the highest confirmed Cooper-pair formation temperature up to 75~K in monolayer FeSe film~\cite{Peng2014} is indeed intriguing and provides a unique opportunity to gain more insight into the origin of high-temperature superconductivity.

The multiple-orbital nature of FeSe, combined with spin and charge degrees of freedom, results in the observation of many interesting phenomena such as nemacity~\cite{Fernandes2014,Yu2018}, orbital-selective Mott transition~\cite{Yi2013,Yu2013,Herbrych2018}, and orbital ordering~\cite{Yi2017}. There are suggestions that the orbital fluctuation may provide a new channel for realizing superconductivity~\cite{Saito2010,Kontani2010}. Direct determination of the correlation between the orbital nature of the low-energy electronic states and the superconducting gap is crucial to understand the superconductivity mechanism of the Fe-based superconductors.

It has been a debate, similar to the cuprate superconductors, whether there exists an antiferromagnetic Mott insulating parent phase for FeSe superconductors. Chen et~al. reported the existence of the $\beta$-Fe$_4$Se$_5$ phase with $\sqrt{5}\times\sqrt{5}$ Fe-vacancy order and argued this $\beta$-Fe$_4$Se$_5$ phase to be the parent phase of FeSe superconductors~\cite{Chen2014}. A theoretical simulation shows that the ground state of $\sqrt{5}\times\sqrt{5}$ Fe-vacancy ordered $\beta$-Fe$_4$Se$_5$ has a pair-checkerboard antiferromagnetic order, and the Fe 3$d$ orbitals govern the low-energy physics~\cite{Gao2017}. Analogously, recent studies demonstrated unambiguously that the K$_2$Fe$_4$Se$_5$, which exhibits a $\sqrt{5}\times\sqrt{5}$ Fe-vacancy order accompanying with an antiferromagnetic order~\cite{Bao2011,Zhao2012}, is the parent phase of the superconductor K$_2$Fe$_{4+x}$Se$_5$~\cite{Wang2015,Wang2019}. Detailed studies of the Fe vacancy in K$_2$Fe$_{4+x}$Se$_5$ reveal that its order/disorder is closely associated with superconductivity.

Fe-deficient $\beta$-Fe$_4$Se$_5$ is expected to exhibit mixed-valence, which may show correlations among the orders of charge, spin, and orbital, resembling magnetite~\cite{Jeng2004,Huang2006}, which exhibits the Verwey transition around 120~K. The metal-insulator transition and the drop in magnetization are the signatures of the Verwey transition~\cite{Walz2002}. In crystal-like $\beta$-Fe$_4$Se$_5$ sheets, a metal-insulator transition and a drop in magnetic susceptibility have been observed at low temperature~\cite{Chen2014}. A more detailed investigation of the Fe-deficient $\beta$-Fe$_{4+\delta}$Se$_5$ will provide clear picture of the physics underlying and contribute to better understanding the superconductivity in FeSe. Nonetheless, in spite of the success of synthesizing FeSe mesoscale materials having various Fe-vacancy orders~\cite{Chang2012,Chang2014,Chen2014}, Fe-vacancy-ordered $\beta$-Fe$_4$Se$_5$ is yet lack of sufficient investigation.

In this work, we present the electrical transport of the Fe-deficient $\beta$-Fe$_{4+\delta}$Se$_5$ single-crystal nanowires with $\sqrt{5}\times\sqrt{5}$-ordered Fe vacancies and mixed valence of Fe. A first-order metal-insulator (MI) transition with an onset transition temperature at $T$~$\sim$~28~K at zero magnetic field is observed. We investigate the magnetotransport and the dielectric relaxation to understand the carrier transport character through the transition. Our study demonstrates the possible topological transition in $\sqrt{5}\times\sqrt{5}$-Fe-vacancy-ordered $\beta$-Fe$_{4+\delta}$Se$_5$.

\begin{figure}[b]
\includegraphics{Fig1
\caption{(a) The SEM image of a Fe$_{4+\delta}$Se$_{5}$ nanowire after deposited electrical leads. The inset is the tilted SEM image, showing that the nanowire has width of $\sim$180~nm and thickness of $\sim$60~nm. (b) The SAED of a Fe$_{4+\delta}$Se$_{5}$ nanowire along the $c$-axis direction. (c) The EELS spectrum of Fe L$_{2,3}$-edges of a Fe$_{4+\delta}$Se$_{5}$ nanowire, showing the mixed valence of Fe.
}
\end{figure}

The $\beta$-Fe$_{4+\delta}$Se$_5$ nanowires were grown by the on-film formation method~\cite{Chang2014}. The nanowires are typically rectangular strips with few-hundred-nanometer width and several-ten-micrometer in length as shown by the scanning electron microscope (SEM) image in Fig.~1(a). We used the 30-nm-thick Si$_3$N$_4$ membranes to hold the nanowires for the benefit of characterizing the nanowires \textit{in situ} by the transmission electron microscope (TEM). The electrical leads, patterned by electron beam lithography with 90~$^\circ\text{C}$ baked PMMA, were made of Au/Ti bilayer (150/30~nm) after the surface oxide layer of the nanowires removed by Ar plasma. The resistances of the nanowires were measured by the four-probe method and lock-in technique with an alternative-current~(AC) excitation less than 300~nA. Sample temperature and external magnetic field are controlled by PPMS (physical property measurement system), and the field orientation to the sample is changed by rotating the sample via a rotation insert.

For sample characterization, we used TEM to obtain the selected area electron diffraction (SAED) pattern, the scanning TEM energy dispersive X-ray spectrometer (STEM-EDS) for element analysis, and the electron energy loss spectroscopy (EELS) to determine the Fe valence. The element mappings show uniform distributions of Fe and Se on the nanowire as shown in Fig.~S1. The ratio between Fe and Se in this nanowire is 4.48/5 according to EDS analysis. As shown in Fig.~1(b),The SAED pattern reveals that the Fe$_{4+\delta}$Se$_5$ nanowire is a single crystal with a tetragonal structure, where the $b$-axis is parallel to the longitudinal axis of the nanowire, the $c$-axis is perpendicular and the $a$-axis is parallel to the wide side of the nanowire. The electron diffraction also clearly shows the super-structure related to the $\sqrt{5}\times\sqrt{5}$ Fe-vacancy order~\cite{Chen2014}. The $a$-axis lattice constant is estimated to be 3.73~{\AA}. The EELS data in Fig.~1(c) demonstrate the absorption peaks for the Fe L$_3$ edge (Fe$^{2+}$ at 706.4~eV and Fe$^{3+}$ at 710.2~eV) and the Fe L$_2$ edge (Fe$^{2+}$ at 718.4~eV and Fe$^{3+}$ at 723.5~eV), confirming the existence of mix-valence state of Fe$^{2+}$ and Fe$^{3+}$.

\begin{figure}[b]
\includegraphics{Fig2
\caption{$B$-field dependence of $R(T)$ with the field direction along the $c$-axis, including results measured in warming and cooling. The upper inset shows $R$ of warming at $T$~$<$~25~K. (b) $T$-$B$ phase diagram includes $T_\text{MI}$ and $T_\text{o}$ derived from $R(T)$ results.
}
\end{figure}

We measured the temperature $(T)$ dependence of the resistance $(R)$ with the temperature control rate of 1~K/min and the AC excitation frequency at 37~Hz. Figure~2 shows the $R$ versus $T$ in cooling and warming measured at different magnetic fields $B$ parallel to the $c$-axis. The transition shows clear thermal hysteresis. At $B$~$=$~0~T, $R$ undergoes a MI transition with the transition onset temperature $T_\text{MI}$~$\sim$~28~K, where $T_\text{MI}$ is defined as the temperature at which the thermal hysteresis appears. As temperature decreases, $R$ drastically increases from $\sim$10~k$\Omega$ and reaches a resistance maximum $\sim$0.6~M$\Omega$ at $\sim$21~K, and then decreases to about 0.4~M$\Omega$ at $T_\text{o}$, which is defined as the temperature where the thermal hysteresis disappears. The thermal hysteresis of $R(T)$ indicates the transition is first-order, which is further confirmed by the applied magnetic field $B$, which enhances the thermal hysteresis loop. The upper inset in Fig.~2 shows more details of the $R(T)$ at $T$~$<$~24~K. A dip is found at $T$~$\sim$~19~K at $B$~$=$~0~T and shifts to higher temperature with increasing $B$ field. The resistance peak value does not change monotonically with $B$ fields. Notice also that the $R(T)$ plateau at $T$~$<$~17~K is solid under $B$ fields up to $B$~$=$~9~T, and reaches a fixed value of $\sim$0.436~M$\Omega$ at 2~K. The lower inset of Fig.~2 displays the phase diagram of the transition between the high-$T$ ($T$~$>$~$T_\text{MI}$) and low-$T$ ($T$~$<$~$T_\text{o}$) phases. Both $T_\text{MI}$ and $T_\text{o}$ increase with increasing $B$.

\begin{figure}[b]
\includegraphics{Fig3
\caption{The AC-excitation frequency, $f$, dependence on the first-order transition. (a) $R(T)$ curves measured in warming at different $f$. The curves of 17~Hz and 37~Hz are almost overlapped. The inset shows the $f$ dependence of conductance $1/R(10\text{K})$ and $1/R(40\text{K})$, normalized by the values at $f$~$=$~17~Hz. The solid lines show $1/R$~$\propto$~$f^{\alpha}$, where $\alpha$ is 1.06 for $1/R(40\text{K})$ and 1.02 for $1/R(10\text{K})$. (b) $1/T_\text{X}(f)$~$-$~$1/T_\text{X}(17\text{Hz})$ for $T_\text{X}$ is $T_\text{MI}$ or $T_\text{o}$. $T_\text{MI}(17\text{Hz})$~$=$~27.4~K and $T_\text{o}(17\text{Hz})$~$=$~17.5~K. The solid lines depict the Arrhenius relation $f(T_\text{X})$~$\propto$~exp$(-E_\text{a}/k_\text{B}T_\text{X})$, where $E_\text{a}$ is 32.3~meV for $T_\text{MI}$ and 44.0~meV for $T_\text{o}$.
}
\end{figure}

To gain more insight into the nature of the transition, we investigated the AC-excitation frequency $(f)$ dependence of $R(T)$ at zero $B$ field. The data were measured from a different section of the same nanowire for the $B$-field dependent experiment. The $R(T)$ at different $f$ are shown in Fig.~3(a). $R$ decreases with increasing $f$ while being weakly $f$-dependent at $T$~$\sim$~28~K. As shown in the inset of Fig.~3(a), as $f$~$>$~37~Hz, $1/R$ shows power-law $f$ dependence, i.e. $1/R$~$\propto$~$f^{\alpha}$ with $\alpha$~$\sim$~1, for both $T$~$<$~$T_\text{o}$ and $T$~$>$~$T_\text{MI}$. The power-law spectrum is known as a fingerprint of dielectric relaxation for hopping systems~\cite{Jonscher1999,Zvyagin2006}. The distinct weak $f$ dependence at $T$~$\sim$~28~K, which is around $T_\text{MI}$ at zero field, suggesting that it is the critical point of the phase transition. The kinks occur in the transition regime at $f$~$\leq$~100~Hz may result from the discrete jumps in order parameters, which has also been observed in the MI transition of mesoscale strongly correlated systems~\cite{Zhai2006,Sharoni2008,Uhlir2016}.

Figure~3(b) shows the $f$ dependence of $T_\text{MI}$ and $T_\text{o}$. When $f$~$>$~37~Hz, $T_\text{MI}$ and $T_\text{o}$ shift with the Arrhenius law $f(T_\text{X})$~$\propto$~exp$(-E_\text{a}/k_\text{B}T_\text{X})$, where $E_\text{a}$ is the activation energy and $T_\text{X}$ is $T_\text{MI}$ or $T_\text{o}$. The Arrhenius law indicates the shift is thermally activated. The extracted $E_\text{a}$'s are 32.3~meV and 44.0~meV for $T_\text{MI}$ and $T_\text{o}$, respectively. The $E_\text{a}$ enhancement across the transition with decreasing $T$ is estimated to be ${\Delta}E$~$=$~$E_\text{a}(T_\text{o})-E_\text{a}(T_\text{MI})$~$\sim$~12~meV. The ratio of $T_\text{MI}$ to ${\Delta}E$ is $\sim$30K/12meV~$=$~$\sim$2.5~K/meV, which is comparable to $\sim$125K/50meV~$=$~$\sim$2.5~K/meV, the ratio of $T_\text{V}$ to the energy gap change deduced from optical conductivity~\cite{Pimenov2005} and photoemission spectrum~\cite{Schrupp2005,Park1997} for single-crystal magnetite. In magnetite studies, the energy gap expansion refers to the increase of $E_\text{a}$ when lowering $T$ over the Verwey transition~\cite{Park1997,Pimenov2005,Schrupp2005,Ramos2008}, and the $E_\text{a}$ values extracted from different experimental approaches~\cite{Kuipers1976,Park1997,Ramos2008,Gooth2014}, including the $f$-dependent shift in AC conductance, are comparable. Thus, ${\Delta}E$ observed in the $\beta$-Fe$_{4+\delta}$Se$_5$ nanowire could be attributed to an energy gap expansion alike to the energy gap expansion observed in the Verwey transition, intimating that charge-orbital ordering may play an important role in the MI transition~\cite{Huang2006,Jeng2004}.

\begin{figure}[b]
\includegraphics{Fig4
\caption{$R$ versus $B$ with $B$ fields parallel to the $c$-axis at different $T$ from 23~K to 35~K (a) and from 17~K to 22~K (b).
}
\end{figure}

As mentioned earlier, the MI transition can be substantially enhanced by applying magnetic fields, as shown in Fig.~2. Consequently, we decided to carry out detailed magnetotransport measurement. Figures~4(a) and 4(b) show the magnetoresistance $R(B)$ of the nanowire at different temperatures with the $B$-field direction parallel to the $c$-axis. At $T$ well above zero-field transition temperature $T_\text{MI}(0\text{T})$, i.e. $T$~$>$~30~K, $R(B)$ exhibits parabolic-like dependence. At $T$ above but close to $T_\text{MI}(0\text{T})$, e.g. $T$~$=$~30~K, $R(B)$ shows linear $B$-field dependence. As 20~K~$<$~$T$~$<$~26~K, a resistance peak in $R(B)$ emerges owing to the $B$-dependent shift of the $R(T)$ peak, and the position of the peak maximum decreases with lowering $T$. Colossal positive magnetoresistance with ${\Delta}R(B)/R(0\text{T})$~$>$~10, where ${\Delta}R(B)$~$=$~$R(B)-R(0\text{T})$, is observed in 26~K~$<$~$T$~$<$~34~K. The colossal magnetoresistance is due to the rise of $T_\text{MI}$ by applying $B$ fields. Strikingly, $R(B)$ at $T$~$\leq$~17~K is essentially independent on $B$. We do not find obvious magnetic hysteresis in magnetoresistance, being consistent with the picture of antiferromagnetism in $\sqrt{5}\times\sqrt{5}$ Fe-vacancy-ordered $\beta$-Fe$_4$Se$_5$~\cite{Chen2014,Fang2016,Gao2017}.

\begin{figure}[b]
\includegraphics{Fig5
\caption{$B$-field orientation dependence. (a) and (b) $R(T)$ measured in warming at $B$ fields of different $\phi_a$ and $\phi_b$, respectively.The insets show ${\Delta}R(B)/R(0\text{T})$ of $T$~$=$~40~K at different $\phi_a$ and $\phi_b$, respectively. (c) and (d) $T_\text{MI}$ and $T_\text{o}$ versus $B_{c}$, the $B$-field component in the direction along the $c$-axis, in $\phi_a$ and $\phi_b$ rotations, respectively. The solid lines represent $T_\text{MI}$ and $T_{o}$ deduced when the $B$ fields parallel to the $c$-axis.
}
\end{figure}

We further investigated the $B$-field orientation dependence of the transition at $B$~$=$~9~T with the nanowire rotated by an angle $\phi_a$ ($\phi_b$) with the rotation axis along the $a$-axis ($b$-axis). The $B$-field direction is parallel to the $c$-axis as $\phi_a$ ($\phi_b$) is zero. The data of $\phi_a$ and $\phi_b$ dependences were acquired from two different sections of the same nanowire. Figures~5(a) and 5(b) show the $R(T)$ results at $B$~$=$~9~T at different $\phi_a$ and $\phi_b$, respectively. As the $B$-field direction changes from parallel to perpendicular to the $c$-axis, the positions of the MI transition and the resistance peak shift to lower $T$. Figures~5(c) and 5(d) show the $T_\text{MI}$ ($T_\text{o}$) versus $B_c$, the $B$-field components parallel to the $c$-axis, where $B_c$~$=$~$9\text{T}\times|\text{cos}(\phi_a)|$ and $B_c$~$=$~$9\text{T}\times|\text{cos}(\phi_a)|$ for different $\phi_a$ and $\phi_b$ at $B$~$=$~9~T, respectively. The $T_\text{MI}$ ($T_\text{o}$) values of different $\phi_a$ and $\phi_b$ scatter along the trace of $T_\text{MI}$ ($T_\text{o}$) versus $B_c$~$=$~$B$, where the $B$-field direction is parallel to the $c$-axis. In other words, $T_\text{MI}$ ($T_\text{o}$) is enhanced anisotropically in $B$ fields with the preferred orientation along the $c$-axis. This observation indicates that the transition may be associated with the spin-orbital coupling in the Fe-vacancy-ordered $\beta$-Fe$_{4+\delta}$Se$_5$ .

The insets in Fig.~5(a) and 5(b) show ${\Delta}R(B)/R(0\text{T})$ at $T$~$=$~40~K of different $B$-field orientations as an example to demonstrate the magnetotransport in the high-$T$ phase. When $B$ fields orient along the $c$-axis, ${\Delta}R(B)/R(0\text{T})$~$\propto$~$B^{\gamma}$, where $\gamma$~$=$~2.7~($\gamma$~$=$~2.6) in the inset of Fig.~5(a)~(5(b)), which is not a conventional Lorentzian magnetoresistance. ${\Delta}R(B)/R(0\text{T})$ decreases when $B$ fields orient from parallel to perpendicular to the $c$-axis. When $B$ fields are perpendicular to the $c$-axis, $\left|{\Delta}R(B)/R(0\text{T})\right|$~$<$~0.01 and is weakly $B$-dependent, showing a two-dimensional-like (2D-like) magnetotransport character and indicating the stripe-shaped nanowire as a quasi-2D system. At $T$~$<$~17~K, $R$ values under applying $B$ fields is rigid regardless of changing $B$-field orientation and magnitude, behaving like one-dimensional (1D) conduction of fully spin polarization~\cite{Chen2008,Chen2012}, reminiscent of a Wigner-like phase~\cite{Wu2012,Zhang2013}. The $B$-field orientation dependence of $R$ peak value is complicated, as shown in Fig.~5(a) and 5(b). When $B$ fields are parallel to the $a$-axis ($b$-axis), the peak height increases (decreases) with increasing $B$. The $R$ peak emerging in the transition is in the phase coexistence regime of the first-order transition. The peak-like feature is also observed in three-dimensional (3D) topological insulators and is regarded as a sign of the existence of the topological surface state~\cite{Cai2018,Wolgast2013,Luo2015,Dankert2018}. The high-$T$ 2D-like phase and the low-$T$ 1D-like phase may be different in topology. $R$ peak feature of the nanowire at $T$~$<$~$T_\text{MI}$ may reveal the emergence of the topological 1D channels. The magnetoresistive transition suggests a possible topological phase transition.

We studied the electrical transport of tetragonal $\beta$-Fe$_{4+\delta}$Se$_5$ nanowires possessing $\sqrt{5}\times\sqrt{5}$ Fe-vacancy ordering and mixed-valance state of Fe$^{2+}$ and Fe$^{3+}$. The first-order MI transition of $T_\text{MI}$~$\sim$~28~K observed at zero $B$ field, and a Verwey-like energy gap expansion of $\sim$12~meV, suggesting that charge-orbital ordering may involve the transition, is extracted from the thermal activation of the dielectric relaxation. When applying external $B$ fields, the MI transition shifts toward higher temperatures remarkably, resulting in the colossal positive magnetoresistance, and the thermal hysteresis gets more pronounced. The $B$-dependent rise of $T_\text{MI}$ and $T_\text{o}$ shows anisotropy with the preferred orientation parallel to the $c$-axis, indicating that the spin-orbital coupling in the nanowire is significant. Based on the results of magnetoresistance, the high-$T$ phase shows 2D-like conduction, and the low-$T$ phase shows the rigidity under external $B$ fields, behaving like fully spin-polarized 1D conduction. The high-$T$ and low-$T$ phases might be topologically distinct, demonstrating the possible topological phase transition in $\sqrt{5}\times\sqrt{5}$ Fe-vacancy-ordered $\beta$-Fe$_{4+\delta}$Se$_5$ nanowires. Our results could provide critical information for better understanding the origin of superconductivity in FeSe-based superconductors and the orbital-related physics in the FeSe-based materials.

\begin{acknowledgments}
We thank Dr.~Peramaiyan~Ganesan and Dr.~Yen-Fu~Liu of the Institute of Physics, Academia Sinica, for valuable discussion. The work is supported by the Ministry of Science and Technology under Grant No.~MOST108-2633-M-001-001 and Academia Sinica Thematic Research Grant No.~AS-TP-106-M01.
\end{acknowledgments}


text:
\section{Introduction}
With the new generation of Imaging Atmospheric Cerenkov Telescopes (IACTs)
like HESS \cite{Aharonian:2006pe}, MAGIC \cite{Albert:2007xh}
and VERITAS \cite{Acciari:2008ah}
Very High Energy (VHE) $\gamma$-ray astronomy
became a very dynamic and exciting discipline with many new detections of
galactic and extra-galactic objects every year.
Almost one third of the objects detected in VHE $\gamma$-rays
are blazars, i.e. Active Galactic Nuclei (AGNs) which contain relativistic jets
pointing approximately in the direction of the observer. Their energy spectra
show no or very weak emission lines, but a continuous distribution with two
broad peaks: one in the UV to soft X-ray band and a second one in the GeV-TeV
range. One of the most interesting aspects of blazars is their variability,
observed in all frequencies and on different time scales ranging from weeks to
minutes \cite{Aharonian:2007ig} \cite{Albert:2007zd}.

The data collected so far in many multiwavelenght observations
is not yet enough to fully constrain the theoretical
models which try to explain the acceleration and emission mechanisms in
blazars. In particular it is not yet clear if the leptonic or hadronic
processes play a decisive role. For example, the Synchrotron-Self Compton
(SSC) leptonic models \cite{Tavecchio:1998xw} can successfully describe most
of the existing data and offer a reasonable explanation of the fast
variability of blazars.
On the other hand, the hadronic models, like the
Synchrotron Mirror Model (SMM) \cite{Reimer:2005sj} or Synchrotron Proton
Blazar (SPB) \cite{Muecke:2002bi} model can also very well describe the SED
structure. Additionally they explain the ''orphan'' gamma-ray flares (like
the one observed for 1ES 1959+650 in 2002 \cite{Krawczynski:2003fq}) and
predict emission of high energy neutrinos.

\section{AGN monitoring}
A long term monitoring of the VHE $\gamma$-ray flux variability performed
in the GeV-TeV energy range will give vital input to understand the problems
mentioned above.

Observations scheduled independently of any knowledge of the
source state provide us with an unbiased distribution of the flux states.
Many of the previously performed measurements were triggered by an observed
enhancement of flux state in other wavebands and therefore
observations of AGNs during a low flux state are still sparse. Consequently
any statistical study which requires high statistics on various flux levels
is difficult. An example of such a study might be the determination of flaring
state probabilities, essential for estimation of
the statistical significance of possible correlations between flaring states
and other observables, such as neutrino events \cite{Satalecka:2007}. Here a long term
monitoring plays a crucial role, especially in a view of the results expected
from the IceCube neutrino observatory \cite{Ahrens:2003ix}.

Another interesting application is the investigation of spectral changes during
periods of different source activity, which can improve our knowledge about the
acceleration and emission processes.

Last but not least AGN monitoring can serve to trigger Target of Opportunity
(ToO) observations. The ToO observations may be performed by the IACT issuing
the ToO trigger but may also include other IACTs - allowing to increase
the time coverage of the observations - or telescopes
and satellites observing other wavelengths. X-ray obsrvations are especially
interesting in a context of the ''orphan'' TeV flares.

\section{The MAGIC Telescope}

MAGIC is currently the largest single-dish IACT
for VHE $\gamma$-ray astronomy. It is located on
the Canary Island of La Palma, at an altitude of 2200\,m a.s.l and has been
in scientific operation since summer 2004.
The camera is equipped with 576 photo-multiplier
tubes (PMTs). In January 2007 a major upgrade of the MAGIC Telescope with
a Multiplexed Fiber-Optic 2\,GSamples/s FADC Data Acquisition system took place
\cite{Goebel:2007cb}. The fast readout minimizes the influence of the background
from the light of the night sky and the additional information on the time
structure of the event helps to reduce the background from hadronic
showers \cite{Aliu:2008pd}. The trigger threshold of MAGIC is around 60\,GeV\footnote{The trigger threshold is defined as the peak of the energy distribution of the triggered events.}.
A source emitting $\gamma$-rays
at a flux level of 1.6\% of the Crab Nebula can be detected with 5\,$\sigma$
significance within 50 hours of observation time. The sensitivity is sufficient
to establish a flux level of about 30\% of the Crab flux above 300\,GeV for a
30\,min observation (Fig.1). A quick on-line analysis system allows to estimate
the flux level of the observed source and facilitate the issuing of ToO triggers
during data taking.
MAGIC can also observe under moderate moon and twilight conditions with only
slightly lower sensitivity.
The construction of a second telescope (MAGIC-II) is being finalized
and the new telescope is undergoing various tests. It is planned
to start stereoscopic observations in spring 2009.
\begin{figure}
  \includegraphics[height=.4\textheight]{AGN_monitoring_fig1.eps}
  \caption{MAGIC sensitivity \cite{LiMa:1983} as a function of exposure time.}
\end{figure}

\section{Monitoring strategy}
Due to the low duty cycle (10-12\%) the operation time of IACTs is very
precious and the observations have to be scheduled very carefully in order to
achieve the right balance between the high priority, deep observations and AGN
monitoring.
Therefore usually previous generation IACTs which are still operational are
used for this purpose. For example the Whipple\footnote{http://veritas.sao.arizona.edu/content/blogsection/6/40/}
telescope, which is leading a monitoring program since 2005 for five well
established TeV sources, or DWARF, one of the former HEGRA telescopes
currently being refurbished and brought back into operation on La Palma \cite{DWARF}. In the case of MAGIC, a high sensitivity latest generation
IACT, some of the monitoring observations can be performed
during moderate moonlight or twilight, which keeps the impact on the overall
observation schedule low.

In order to achieve a dense sampling, up to 40 short observations per source
are scheduled, evenly distributed over the observable time during MAGIC
observation period. Each of them should be long enough to detect a given minimum
flux level taking into account the sensitivity of the telescope. In the case
of MAGIC, three sources were chosen for a regular monitoring: Mrk 421, Mrk 501
and 1ES 1959+650. The first two are relatively bright and usually 15-30\,min
observations are scheduled for them. 1ES~1959+650 being fainter requires
longer observation times, at least 30 minutes per single exposure.

\section{Results}
In this section we will present the preliminary results of the MAGIC AGN
monitoring program for the observation season 2007/2008. The data have been
processed with the standard MAGIC analysis tools \cite{Albert:2007xh}.
A fraction of the data has been removed due to poor observation conditions.
All cuts are optimized and verified with Crab data
in the corresponding observational seasons.
\subsection{Mrk 421}
Mrk 421 is the first extragalactic object detected at energies above 500\,GeV
\cite{Punch:1992} \cite{Petry:1996yj} and one of the best studied.
Here the preliminary results of the data which were
taken from February 2007 to June 2008 are presented (82 hours).
After the data selection, about 66 hours of data
(80\%) were used for further analysis. The observations were mostly
performed in wobble mode, which allows to simultaneously collect signal
and background events.

Fig.\ref{LC_Mrk421_all} shows all the data collected from Mrk 421 by
the MAGIC telescope since 2004 \cite{Albert:2006jd} \cite{Goebel:2007uu}.
In Fig.\ref{LC_Mrk421_0708} the 2007/2008 observational period is shown magnified.
We can clearly see that Mrk 421 was especially active in 2008 and
the flux was never below 0.5\,Crab. During this time of high activity,
in accordance
with the ToO agreement, MAGIC issued several alerts to HESS and VERITAS.
Detailed analysis of the collected data, such as intra-night variability and
spectrum measurement will be published elsewhere.

\begin{figure}
  \includegraphics[width=1.\textwidth]{AGN_monitoring_fig2.eps}
  \caption{\label{LC_Mrk421_all}
Mrk 421 light curve showing all data collected by MAGIC \cite{Albert:2006jd} \cite{Goebel:2007uu}.
}
\end{figure}

\begin{figure}
  \includegraphics[width=1.\textwidth]{AGN_monitoring_fig3.eps}
  \caption{\label{LC_Mrk421_0708}
Mrk 421 light curve showing MAGIC data from AGN monitoring season 2007/2008.
}
\end{figure}

In the SSC framework, a correlation between TeV photon and the other
wavelengths, especially X-ray, is expected \cite{Maraschi:1999pg}.
Since usually the TeV
data is relatively poorly sampled, it is not easy to get really simultaneous
data, we therefore searched for X-ray and optical measurements which are
within 6\,h before or after the TeV observations.

In this work we use the X-ray data taken by the RXTE/ASM satellite, which is
designed for all sky monitoring observations. The data are available
on the web page\footnote{http://xte.mit.edu/ASM\_lc.html} immediately after the
observations. The averages and errors of ASM data points are calculated
on a dwell-by-dwell\footnote{one dwell lasts 90 seconds} data basis.
If the number of dwells are fewer than five,
we discard that data point; 75 TeV-X-ray measurements pairs were
finally selected (shown in Fig.\ref{TeVXray}).

The optical R-band data were provided by the Tuorla Observatory
Blazar Monitoring Program\footnote{http://users.utu.fi/kani/1m/} with the
1.03\,m telescope at the Tuorla Observatory Finland and the 35\,cm KVA telescope
at La Palma, Canary Islands. From the optical flux the flux of the host galaxy
and the flux contribution from the companion galaxies were substracted
\cite{Nilsson:2007ax}. Eventually, 56 measurements pairs were found
(Fig.\ref{TeVoptical}).

The correlation coefficients for both plots were calculated using Pearson's method
(see section 14.5 in \cite{638765}).
For X-ray/TeV correlation, the correlation coefficient value is r~ = 0.77$\pm$0.05.
For optical/TeV plot, the coefficient value is r = 0.03$\pm$0.14.
We thus found a significant correlation with X-rays.
This type of correlation is easier accommodated in leptonic models
(see e.g. \cite{Wagner:2008cw} and references therein) but there also exist
certain hadronic models which support such correlation (see e.g.
\cite{Muecke:2002bi} or \cite{Aharonian:2000pv}).
However, one has to keep in mind, that the data were not taken strictly
simultaneously and any conclusion as to the origin or mechanisms of the emission
should be made very carefuly.

\begin{figure}
  \includegraphics[height=.5\textheight]{AGN_monitoring_fig4.eps}
  \caption{\label{TeVXray}VHE $\gamma$-ray (MAGIC) - X-ray (ASM) correlation plot for Mrk 421.}
\end{figure}

\begin{figure}
  \includegraphics[height=.5\textheight]{AGN_monitoring_fig5.eps}
  \caption{\label{TeVoptical}VHE $\gamma$-ray (MAGIC) - R-band (Tuorla) correlation plot for Mrk 421.}
\end{figure}

\subsection{Mrk 501}
The BL Lac object Mrk 501 was the second established TeV-blazar
\cite{Quinn:1996dj} \cite{Bradbury:1996mi}. MAGIC has been observing this source
since the year 2005 \cite{Albert:2007zd}. Here we present preliminary results
based on an observation time of 16 hours (after quality selection), collected
between February 2007 and August 2008. In order to maximize the time coverage
of this source, observations were carried out mostly in the presence of moderate
moonlight or twilight (9 hours, i.e. 56\% of the observation time).
The observations were performed in wobble mode.

Fig. \ref{LC_Mrk501_all} shows all the data collected form Mrk 501 by the MAGIC
telescope since 2005  \cite{Albert:2007zd} \cite{Goebel:2007uu}.
In Fig. \ref{LC_Mrk501_0708} the 2007/2008 monitoring season is shown magnified.
During this period the source was in a relatively low state
(below 1\,Crab). On the other hand a dense sampling was obtained and a
statistical analysis, e.g. estimation of the of the source state probability,
based on the above data is in progress.

\begin{figure}
  \includegraphics[width=1.\textwidth]{AGN_monitoring_fig6.eps}
  \caption{\label{LC_Mrk501_all}
Mrk 501 light curve showing all data collected by MAGIC \cite{Goebel:2007uu} \cite{Albert:2007zd}.
}
\end{figure}

\begin{figure}
  \includegraphics[width=1.\textwidth]{AGN_monitoring_fig7.eps}
  \caption{\label{LC_Mrk501_0708}
Mrk 501 light curve showing MAGIC data from AGN monitoring season 2007/2008.
}
\end{figure}

\subsection{1ES 1959+650}
The observation of an ''orphan'' TeV flare (i.e. not accompanied by a
corresponding X-ray
flare) in 2002 \cite{Krawczynski:2003fq} made 1ES 1959+650 one of
the most interesting VHE $\gamma$-ray sources.
MAGIC has been observing this object
since 2004 \cite{Aliu:2005jv}, a regular monitoring was first performed
in the 2005/2006 season \cite{Goebel:2007uu}.
Unfortunately, the source position allows observations only under a relatively
high zenith angles (35-50\,deg).
Moreover, when not in a high state the source is rather faint compared to
those discussed above.

MAGIC observed 1ES 1959+650 from April 2007 till October 2008 for 27 hours
(including 13 hours of data taken during moonlight or twilight).
Observations were performed in wobble mode. After rejection of data due to poor
quality 3.4 hours of data taken in 2008 were analyzed. The overall significance
of this sample is 2.5\,$\sigma$, which allows us to set an upper limit on the flux
above 300\,GeV of 2.53\,$\times 10^{-11}$\,ph/cm$^2$/s at 90\% C.L. Comparing to
previous results \cite{Aharonian:2003} \cite{Gutierrez:2006ak} we can state
that MAGIC observed this source in its usual low state.
\section{Conclusions}
During the observational season 2007/2008 MAGIC monitored three sources:
Mrk 421, Mrk 501 and 1ES 1959+650. The detailed analysis of the collected data
is going on. In this work we presented preliminary results of the measured
flux levels for all three sources. We also investigated the possible
correlation of the TeV and X-ray/optical flux levels for Mrk 421. We found a
significant correlation with X-rays and no correlation with optical R-band.
The most interesting material for further studies are the
many fast flares observed for Mrk 421 in 2008. Mrk 501 and 1ES 1959+650 were
found in low state, but dense sampling of Mrk 501 flux states looks promising
for statistical studies.

\begin{theacknowledgments}
  We would like to thank the Instituto de Astrofisica de Canarias for the
excellent working conditions at the Observatorio del Roque de los Muchachos
in La Palma. The support of the German BMBF and MPG, the Italian INFN and
Spanish MCINN is gratefully acknowledged. This work was also supported by
ETH Research Grant TH 34/043, by the Polish MNiSzW Grant N N203 390834,
and by the YIP of the Helmholtz Gemeinschaft.
\end{theacknowledgments}

\bibliographystyle{aipproc}