\section{The RF-dressed shell trap}
\label{sec:ShellTrap}

We consider the case of a $^{87}$Rb atom in the presence of a static quadrupole field
$\boldsymbol{B}_\text{q}\left(\mathbf{r}\right)=\alpha (x\hat{\mathbf{e}}_{\text{x}}+y \hat{\mathbf{e}}_{\text{y}}-2 z\hat{\mathbf{e}}_\text{z})$, where $\alpha$ is a magnetic quadrupole gradient and $(\hat{\mathbf{e}}_{\text{x}},\hat{\mathbf{e}}_{\text{y}}, \hat{\mathbf{e}}_\text{z})$ are the unit vectors. A homogeneous RF-field of arbitrary polarization $\boldsymbol{B}_{\textrm{RF}}(t)$ and frequency $\omega_\text{RF}$ is also applied. The Hamiltonian for weak fields $|g_\text{F} \mu_\text{B} \boldsymbol{B}_\text{q}(\mathbf{r})|\ll \hbar \omega_\text{hfs}$ (where $\omega_\text{hfs}$ is the hyperfine splitting frequency, and $\mu_{B}$ is the Bohr magneton) reads:
\begin{equation}
H = \frac{A_{\text{hfs}}}{\hbar^{2}} \boldsymbol{I}\cdot\boldsymbol{J} + \frac{\mu_\text{B}}{\hbar}(g_I \boldsymbol{I} + g_J \boldsymbol{J}) \cdot (\boldsymbol{B}_\text{q}(\mathbf{r})+\boldsymbol{B}_{\textrm{RF}}(t)),
\label{eq:H_uncoupled}
\end{equation}

\noindent where $A_{\text{hfs}}$ is the hyperfine constant of $^{87}$Rb, $\boldsymbol{I}$ and $\boldsymbol{J}$ are the nuclear and electronic angular momentum operators, with $g_{I}$ and $g_{J}$ being their respective gyromagnetic ratios. For weak magnetic fields, the total angular momentum  $\boldsymbol{F}=\boldsymbol{I}+\boldsymbol{J}$ with eigenstates in a coupled basis $|F,m_{F}\rangle$ can be used instead of $|I,m_{I},J,m_{J}\rangle$. For $^{87}$Rb at the electronic ground state ($I=3/2$, $J=1/2$) there are two hyperfine subspaces with $F=I\pm 1/2$, and hyperfine Zeeman sub-states given by $m_{F}=-F\dots 0 \dots F$ in each one. The total dimension is $(2J+1)(2I+1)=8$. The total angular momentum gyromagnetic factor $g_{F}$ is (see \cite{Steck2001} and references therein)

\begin{equation}
g_{F}=g_{J}\frac{F(F+1)-I(I+1)+J(J+1)}{2F (F+1)}+g_{I}\frac{F(F+1)+I(I+1)-J(J(J+1))}{2F(F+1)}
\label{eq:gFdefinition}
\end{equation}

\noindent with $g_{J}=2.002331$ and $g_{I}=-0.000995$. Eq.~\ref{eq:gFdefinition} results in different $g_{F}$ for $F=1,2$ ($g_{1}=-0.5018$ and $g_{2}=0.4998$) \cite{Steck2001}.
Eq.~\ref{eq:H_uncoupled} is generally solved and simplified by applying the Rotating Wave Approximation (RWA) and by using a \textit{dressed} basis $|F,\bar{m}_{F}\rangle$, of the same dimension as $|F,m_{F}\rangle$. The result, including gravity, is the adiabatic dressed potential~\cite{Zobay2001}:

\begin{equation}
V^{F}\left(\mathbf{r}\right)=s\left(I+\frac{1}{2}\right)\frac{\hbar \omega_\text{hfs}}{2}+ s\bar{m}_{F} \hbar  \sqrt{\delta_{F}^{2} + \Omega_{F}^{2}\left(\mathbf{r}\right)}+M g z
\label{eq:dressedPotentiala}
\end{equation}

\noindent where $s=g_{F}/|g_{F}|$ is the sign of $g_{F}$, with
  $s=-1$ for $F=1$ and  $s=+1$ for $F=2$.
  $M$ is the atomic mass of $^{87}$Rb and $g$ is the earth's gravitational acceleration.
The detuning of the RF from the resonance is $\delta_{F}=\Omega_\text{L}^{F}\left(\mathbf{r}\right)-\omega_{F}$ where $|\hbar \Omega_\text{L}^{F}(\mathbf{r})|= |g_{F} \mu_{\text{B}} \boldsymbol{B}_{q}(\mathbf{r})|$ is the Larmor frequency and $\omega_{F}$ is the RF-frequency $\omega_\text{RF}$ that dresses each manifold, i.e. \,$\omega_{1}$ for F=1 and $\omega_{2}$ for F=2. The Rabi coupling, $\Omega_{F}(\mathbf{r})$, varies with the position and depends on the spin gyromagnetic ratio \cite{Garraway2015}. Please note, that the sign of $g_{F}$ determines which local polarization component of the RF couples the $m_F$ states, i.e.~for the absorption of one rf-photon $(\Delta m_F=+1)$ one requires $\sigma_{-}$ for the $F=1$ and $\sigma_{+}$ for the $F=2$ state.
Since the $F=1$ and $F=2$ manifolds are interacting with mutually orthogonal RF-polarizations, they can be addressed entirely independently. A circularly polarized RF-field in the laboratory frame can be written as:

\begin{equation}
\boldsymbol{B}_\text{RF}=B_{\text{RF}}\left(\cos{\left(\omega_{F} t\right)}\hat{\mathbf{e}}_{\text{x}}+s\sin{\left(\omega_{F} t\right)}\hat{\mathbf{e}}_{\text{y}}\right)
\label{eq:monochrombfield}
\end{equation}

\noindent Each sense of rotation, specified with $s=\pm 1$ and defined with respect to the symmetry axis of the quadrupole field, couples atoms to either of the F manifolds. The local coupling with respect to the spatially dependent quantization axis defined by the direction of $\boldsymbol{B}_{q}(\mathbf{r})$ is:

\begin{equation}
|\Omega_{F}\left(\mathbf{r}\right)|=\frac{\Omega_{\text{RF}}}{2}\left(1-\frac{2z}{\sqrt{\rho^{2}+4z^{2}}}\right)
\label{rabicoupling}
\end{equation}

\noindent where $\hbar \Omega_{\text{RF}}=|g_{F}|\mu_{\text{B}}B_{\text{RF}}$ is the maximum Rabi frequency (note that all axes have an origin at the zero of the quadrupole field). The resulting potential is an iso-magnetic surface defined by the resonance of $\omega_{F}$ and $\boldsymbol{B}_{q}\left(\mathbf{r}\right)$ \cite{Lesanovsky2006}. Trapping requires low-field-seekers, i.e.~states with $g_{F} m_{F}> 0$, when $\gamma=\frac{M g}{2 m_{F} g_{F} \mu_{\text{B}} \alpha}<1$ \cite{Merloti2013}.
The trapping frequencies close to the trap center are:

\begin{eqnarray}
\omega_{\rho}= \sqrt{\frac{g}{4 z_{0}}}\left(1-\frac{|m_{F}|\hbar \Omega_{\text{RF}}}{2 m g z_{0}}\sqrt{1-\gamma^{2}}\right)^{\frac{1}{2}} \label{eq:shelltrapfreqrad} \\
\omega_\text{z}=\frac{2 |g_{F}| \mu_{B} \alpha}{\hbar} \sqrt{\frac{|m_{F}| \hbar}{M \Omega_{\text{RF}}}}\left(1-\gamma^{2}\right)^{\frac{3}{4}}
\label{eq:shelltrapfreqaxial}
\end{eqnarray}

\noindent where:

\begin{equation}
z_{0}=\frac{\hbar \omega}{2\alpha \mu_{\text{B}} |g_{F}|}\left(1+\frac{\gamma}{\sqrt{1-\gamma^{2}}}\frac{\Omega_{\text{RF}}}{\omega}\right)
\label{eq:trapposition}
\end{equation}
is the vertical position of the trap. The values of $\gamma$ for the the spin states $|1,\bar{m}_{F}=-1\rangle$ and $|2,\bar{m}_{F}=1\rangle$ states differ only by a small amount, e.g. $\gamma_{1}/\gamma_{2}=g_{2}/g_{1}\approx0.996$.
This small difference causes the two traps to have both different positions $z_0$  and trapping frequencies.
Therefore, any superposition between the two states  will  dephase rapidly, thus  severely limiting any interferometric measurement using these two states.

\section{The Bi-chromatic Shell Trap \label{sec:BiChromaticShellTrap}}

As stated above, the two manifolds $F=1$ and  $F=2$ are dressed by  RF-fields, which are mutually orthogonal and can have different frequencies and/or amplitudes.
The  mismatch of the traps can be overcome at least partially using \emph{bi-chromatic} adiabatic  potentials.
The RF-dressing fields used in the bi-chromatic shell trap are composed of one  $\sigma_{-}$ polarized field with $(\omega_{1}, \Omega_{1})$ and one $\sigma_{+}$ polarized field with $(\omega_{2},\Omega_{2})$. They can be generated from  rf-coils placed in the x- and y-direction with the linearly polarized fields as:

\begin{equation}
\begin{split}
\boldsymbol{B}_\text{\text{RF}}(t)=
\boldsymbol{B}_{1}(t)+\boldsymbol{B}_{2}(t)=
B_{\text{x1}} \cos{\left(\omega_{1}t\right)}\hat{\mathbf{e}}_{\text{x}}-
B_{\text{y1}} \sin{\left(\omega_{1}t\right)}\hat{\mathbf{e}}_{\text{y}}
\\+
B_{x2} \cos{\left(\omega_{2}t\right)}\hat{\mathbf{e}}_{\text{x}} +
B_{y2} \sin{\left(\omega_{2}t\right)}\hat{\mathbf{e}}_{\text{y}}
 \end{split}
 \label{eq:bifield}
 \end{equation}

 \noindent which results in two circularly polarised components at the position of the atoms when the amplitude of both frequency components at the RF-coils is the same.
 %

The traps for the two states have independently tunable position and trapping frequencies.
Perfect matching would be achieved when the positions and both trapping frequencies are identical for both states at the same time.
Unfortunately, this cannot be achieved simultaneously  for both states, see Eqs.(\ref{eq:shelltrapfreqrad},\ref{eq:shelltrapfreqaxial},\ref{eq:trapposition}).
It is possible, however, to  perfectly match $z_{0}$ and $\omega_{z}$ and at the same time keep the difference in $\omega_{\rho}$
very small.
For a given $\alpha$ and by fixing the parameters for one trap ($\omega_1$ and $\Omega_1$), we can calculate  the ones for the other:

\begin{eqnarray}
\Omega_{2}=\left(\frac{g_{2}}{g_{1}}\right)^{2}\left(\frac{1-\gamma_2^2}{1-\gamma_1^2}\right)^{3/2}\Omega_1
\label{eq:matchingcondrabi}\\
\omega_{2}=\frac{g_{2}}{g_{1}}\left(\frac{\gamma_{1} \Omega_{1}}{\sqrt{1-\gamma_{1}^2}}+\omega_{1}\right)-\frac{\gamma_2 \Omega_2}{\sqrt{1-\gamma_{2}^2}}\label{eq:matchingcondfreq}
\end{eqnarray}

\noindent The full analytic expression for the resulting radial frequency difference is too complex to be printed here.
Fig.~\ref{fig:radialfreqdiff} shows the difference in radial frequencies $\Delta \omega_{\rho}/2\pi$  for a range of Rabi and RF frequencies assuming the matching conditions of Eq.(\ref{eq:matchingcondrabi}).
It is clear that the frequency difference decreases with  increasing radio frequency $(\omega_{\textrm{RF}})$ and decreasing gradient $(\alpha)$ and Rabi frequency $(\Omega_{\textrm{RF}})$.
We find, for example, that $\Delta \omega_{\rho}/2\pi$ can be as low as $10$\,mHz for Rabi frequencies of the order of $100$\,kHz, RF-frequencies of the order of 2.5MHz and gradients of the order of $20$\,Gpcm.
Assuming only a single vibrational level is occupied, this implies a de-phasing time of the order of 10\,s.

\begin{figure}
\begin{center}
\resizebox{0.8\textwidth}{!}{%
  \includegraphics{gpcm_radial_freq_diff2.pdf}}
  \caption{Difference in radial frequency between the two traps, calculated as $\Delta \omega_{\rho}/2\pi$. For any $B_q$, the position $z_0$ and the axial trapping frequency $\omega_{z}$ are perfectly matched. One can match $\omega_{r}$ but only for asymptotically non-physical parameters, e.g. for very large RF-frequencies and for very small Rabi frequencies and quadrupole gradients.}
  \label{fig:radialfreqdiff}
  \end{center}
\end{figure}

\section{Microwave beam splitters for adiabatic potentials \label{sec:BeamSplitters}}

\iffalse
BS
Fluctuations
 DC
-- Homogeneous
-- QP
AC
-- RF
-- mw
\fi

In the atom interferometer based on the adiabatic traps \cite{Navez2016}, beam splitters are microwave transitions between the clock states.
This adds a term proportional to $\boldsymbol{B}_\text{MW}(t)$ to the total magnetic field term in Eq.~\ref{eq:H_uncoupled}.
The solution of this Hamiltonian for a weak MW-field that acts as a probe on the RF-dressed energy levels leads to inter-manifold transitions \cite{Haroche1970,Nottingham,Luksch2018} on a spectrum of 7 groups (spaced by the RF-dressing frequency) of 5 transitions (spaced by the RF-dressing Rabi frequency), when the initial state is $|F=1,\bar{m}_{F}=-1\rangle$.
These are given by the condition:
$
	\omega_{\text{MW}} = \omega_\text{hfs} + k\,\Omega_{\text{RF}}-n\,\omega_\text{RF} ,
	\label{eq:resonanceCondition}
$
with $n=-3...0...3$ and $k=-3...0...1$.
For all $n$ and $k=0$ transitions between $|1,-1\rangle \rightarrow |2,1\rangle$ occur. Note that we have written this simplified formula for the case where both manifolds are dressed by Rabi frequencies of equal magnitude (see \cite{Nottingham} for a detailed study).
Concretely, the transition with $n=0,k=0$ is coincident with the hyperfine splitting \cite{Nottingham}.
In practice, this transition is shifted both by the non-linear dependence of the eigen-energies on the magnetic field and on the difference in magnitude of $g_1$ and $g_2$.
In addition, we find that the cross-coupling between $B_{2}$ ($B_{1}$) and $F=1$ ($F=2$) induced by the polarization defects results in a periodic modulation of the energy splitting that causes multi-photon transitions \cite{Silveri2017, Oliver2005}, and also time-averaged reshaping of the traps.
In the ideal case, the bi-chromatic shell should exhibit the spectral response explained above, i.e.\,7 peaks close to $\omega_\text{MW}=\omega_\text{hfs}+n\omega_\text{RF}$.
However, we observe a finite number of harmonics separated by $\Delta \omega_{F} =|\omega_{1}-\omega_{2}|$.
We explain these multi-photon transitions with a semi-classical two-level toy model in \ref{app:multiphoton}.

In Section~\ref{sec:BiChromaticShellTrap} we found an approximation for the $\Delta \omega_F$ that matches the shell traps of the two states.
In practice, however, the exact matching is affected by non-linear Zeeman shifts and the experimental  imperfections in the generation of the (RF-)field.
In this section, we present a method to match the two traps experimentally  by spectroscopic means: we choose a pair of frequencies $\omega_{1}$, $\omega_{2}$ close to the prediction of Eqs.~(\ref{eq:matchingcondrabi},\ref{eq:matchingcondfreq}) for $\alpha=45$\,G/cm and $\Omega_{1}/2\pi\approx \Omega_{2}/2\pi \approx 250$\,kHz (see Appendix~\ref{app:traploading} for an explanation of the sample preparation and trap loading).
We then find the longest coherence times by measuring the line-width of the transition  for different values of $\Delta \omega_\text{F}$. To do this, we load the bi-chromatic shell trap for several pairs of $\omega_{1}$ and $\omega_{2}$ and then we fit Lorentzian curves to the number of atoms transferred to the upper state as we scan the microwave frequency $\omega_{\text{MW}}/2 \pi$. By finding the combination of $\omega_{1}$ and $\omega_{2}$ that gives the minimum linewidth, we can thus match the two traps (see Fig.~\ref{fig:tuningfrequency} a)), at which point  we observe a reduction in line-width by of the order of $\times 10$ with respect to the mono-chromatic shell (measured elsewhere \cite{Nottingham}).

In the mono-chromatic shell trap we observe very fast decay curves and no Rabi oscillations, and a line width of $\Delta \nu\approx 1000$\,Hz.
In contrast to this, in the bi-chromatic shell we measured a line-width of $\Delta \nu \approx100$\,Hz for $\Delta t_\text{MW}=10$\,ms, and we observe Rabi oscillations with a decay time of $>10$\,ms. Fig.~\ref{fig:tuningfrequency} b) shows an example of the Rabi oscillations in the bi-chromatic shell trap, with $\alpha =50$\,G/cm, $\omega_{1}/2 \pi=2.294336$\,MHz, $\omega_{2}/2 \pi=2.285238$\,MHz and $\Omega_{1}/2\pi\approx \Omega_{2}/2\pi\approx 270$\,kHz. In the figure, $f_{2}$ is the normalized atom number in the dressed $|2,1\rangle$ state. We fit the data to the the following function $f_{2}=f_{o,2}+(f_{s,2}/2)\left(1- e^{-t/\tau}   \cos{\left(\Omega_{\textrm{MW}} t_{\textrm{mw}}\right)}\right)$. In this equation, $\Omega_{\textrm{MW}}$ is the MW driving field Rabi frequency in the RF-dressed TLS, and $\tau$ is the decay time of the Rabi oscillations amplitude for this transition.
A fit to the data in Fig.~\ref{fig:tuningfrequency} b) yields $\Omega_{\textrm{MW}}/2\pi=146\pm 1$\,Hz and $\tau=18\pm 2$\,ms. The term $f_{0}$ is an offset that accounts for atoms that are previously transferred to the $|2,1\rangle$ due to an imperfect switching on of the MW radiation. The fields at $\omega_{1}$ and $\omega_{2}$ are produced from two separate two-channel sources and then mixed before the amplifiers that are connected to the RF-antennae in the x- and y-direction (see~\ref{app:rfgeneration}).

\begin{figure*}
    %
          \begin{minipage}{.47\textwidth}\centering
     \subfloat[][]{\includegraphics[height=1.9in]{180723_1111HMLndwdth.pdf}\label{fig:tuningfrequencya}}
     \end{minipage}
          \begin{minipage}{.47\textwidth}\centering
     \subfloat[][]{\includegraphics[height=1.9in]{rabioscillations.pdf}\label{fig:tuningfrequencyb}}
     \end{minipage}%
     %
     %
     %
     \caption{(a): measured line-width $\Delta \nu$ in the $n=1$ transition with a microwave pulse of $10$\,ms for the clock transition in the dressed states for different values of $\Delta \Omega$. $\Omega_{1}/2\pi\approx\Omega_{2}/2\pi\approx  250$\,kHz, $\alpha=45$\,G/cm. (b): Example of Rabi oscillations in the bi-chromatic adiabatic potential.}
\label{fig:tuningfrequency}
\end{figure*}

\section{Dephasing and decoherence\label{sec:DephasingSources}}

The RF frequencies are generated with reference to an atomic clock and as such do not contribute to the dephasing.
The energy levels of the atoms in the shell trap are insensitive to changes in the homogeneous DC-fields since their only effect is to move the center of the quadrupole trap.
However, they are  sensitive to fluctuations in the gradient of the quadrupole field and in the Rabi frequencies of the RF-fields, especially to the difference of the Rabi frequencies of the two states, i.e.\,$\Delta \Omega= \Omega_{2}-\Omega_{1}$.
The transition is relatively insensitive to fluctuations in the Rabi frequency as long as $\Delta \Omega$ remains stable.
Fig.~\ref{fig:noisesourcesfig} a) shows the sensitivity of the transitions with respect to fluctuations of the magnetic gradient field.
We write them as $\partial \Delta E/\partial \alpha$, where $\Delta E=V^{2}-V^{1}$ is the difference in energy between the two adiabatic potentials retrieved from Eq.~\ref{eq:dressedPotentiala}, assuming the traps are matched in position and axial trapping frequency at all times following Eqs.~(\ref{eq:matchingcondfreq},\ref{eq:matchingcondrabi}).
Clearly, this sensitivity increases with decreasing gradients.
Fig.~\ref{fig:noisesourcesfig} b) shows the sensitivity of the transition on fluctuations in the difference of Rabi frequencies $\Delta \Omega$, again assuming that the traps are perfectly matched and that the fluctuations in the two RF-fields are uncorrelated.
In this case, the sensitivity is smaller for lower quadrupole gradients $\alpha$.

\begin{figure}
 \includegraphics[width=1 \linewidth]{figures_noise.pdf}
\caption{ (a) We show the sensitivity of the transition to fluctuations in the gradient for various gradients and various $\Omega_{2}$ for fully matched traps (i.e.\,Eqs.~(\ref{eq:matchingcondrabi},\ref{eq:matchingcondfreq}) are fulfilled). (b) We show the susceptibility of the transition associated to fluctuations in the difference in Rabi frequency, that is $\Delta \Omega= \Omega_{2}-\Omega_{1}$.
}
\label{fig:noisesourcesfig}
\end{figure}

\subsection{Sensitivity to homogeneous fields
\label{sec:LineShiftWithBiasField}}

Compared to other magnetically trapped atom-clock interferometers, the shell trap has the considerable advantage that homogeneous DC magnetic field fluctuations do not influence the transition, thus eliminating the need for DC magnetic shielding.
The reason behind this is that the shell trap is based on a magnetic quadrupole field.
Any homogeneous field simply shifts the trap center without affecting its trapping parameters.
In order to check that this shift does not affect the trap, e.g.~due to field inhomogeneities, we measured the clock frequencies of the dressed  $|1,-1\rangle\rightarrow|2,1\rangle$ transition (by finding the center of Lorentzian fits to the atomic population in the dressed $|2,1\rangle$ after short pulses) for extreme variations in homogeneous background field ($0 \leq B_{\text{z}} \leq 3.5$\,G).
For this we loaded a bi-chromatic shell with $\omega_{1}/2\pi=2.294336$\,MHz, $\omega_{2}/2\pi=2.285238$\,MHz, $\Omega_{1}/2\pi\approx \Omega_{2}/2\pi\approx 240$\,kHz and $\alpha=70$\,G/cm (see~\ref{app:traploading}).
Then we drove the transition at $n=0, k=0$ ($\Delta t_{\text{MW}}=2$\,ms) for various field strengths $B_{\text{z}}$.
Fig.~\ref{fig:bfieldshift} shows the detuning $\delta_{\bar{\omega}_\text{MW}}=(\omega_{\text{MW}}-\bar{\omega}_{\text{MW}})/2 \pi$, with the shaded yellow area corresponding to the confidence interval given by the statistical uncertainty of a linear fit to the data points. Here $\bar{\omega}_{\text{MW}}$ is the mean value of all measurements of the transition frequency $\omega_{\text{MW}}/2\pi$ for different homogeneous fields $B_{z}$.
This agrees with the experimental observation that there is no shift in the transition frequency of $|1,\bar{m}_{F}=-1\rangle \rightarrow |2,\bar{m}_{F}=1\rangle$ in the range from zero to 3.5\,G. When we fit a linear slope to the set of data in Fig.~\ref{fig:bfieldshift}, we find a near zero slope ($2\pm10$\,Hz/G) and zero offset ($4\pm24$\,Hz).

\begin{figure}
\begin{center}
\resizebox{0.5\textwidth}{!}{%
  \includegraphics{biasFieldShift.pdf}}
\caption{Measured $\delta_{\bar{\omega}_\text{MW}}=(\omega_{\text{MW}}-\bar{\omega}_{\text{MW}})/2 \pi$ where $\bar{\omega}_\text{MW}/2\pi=\omega_\text{hfs}+2\pi\, 6210$\,Hz
is the mean value of the measured transition frequencies.
The trap parameters were $\omega_{1}/2\pi=2.294336$\,MHz, $\omega_{2}/2\pi=2.285238$\,MHz, $\Omega_{1}\approx \Omega_{2}\approx240$\,kHz, $\alpha=70$\,G/cm. The MW pulse lasted $\Delta t_{\text{MW}}=2$\,ms, much less than one Rabi cycle.
We fit a linear slope to this data set and find it to be $(2\pm10)$\,Hz/G, with an offset of $(-4\pm24)$\,Hz. The yellow area represents the confidence interval given by the statistical uncertainty of the fit.
}
\label{fig:bfieldshift}
\end{center}
\end{figure}

\begin{figure}
\begin{center}
\resizebox{0.5\textwidth}{!}{%
  \includegraphics{181119_1006HM_shiftDeltaOmega.pdf}}
\caption{Measured shift for the dressed $|1,\bar{m}_{F}=-1\rangle \rightarrow |2,\bar{m}_{F}=1\rangle$ transition for several values of the difference in Rabi frequency $\Delta \Omega=\Omega_{2}-\Omega_{1}$ between the $V^{1}$ and $V^{2}$ shell traps.
The solid red line is a fit that yields a slope of $1.01\pm 0.03$\,Hz/Hz with an offset of $-4412\pm 20$\,Hz.
The statistical error in the frequencies is much smaller than the dot-size.
There is a systematic error of about $ 1$\,kHz due to the uncertainty of the determination of $\Delta \Omega/2\pi$ .
}
\label{fig:DeltaOmegaShift}
\end{center}
\end{figure}

\subsection{Line shift with $\Delta \Omega$}
\label{sec:shiftRabi}

We can estimate from Eq.\,\ref{eq:dressedPotentiala} the line shift with the difference in Rabi frequency ($\Delta \Omega$) under the assumption that the position of the two traps is matched.
This shift becomes directly proportional to $\Delta \Omega$ when the gradient is high (see Fig.~\ref{fig:noisesourcesfig} b)).
Fig.\,\ref{fig:DeltaOmegaShift} shows the transition frequency $\delta_\text{MW}=(\omega_\text{MW}-\omega_\text{hfs})/2 \pi$ for several $\Delta \Omega$.
The Rabi frequencies had been previously calibrated by measuring the transitions induced by an additional weak RF-field \cite{GarridoAlzar2006}. To do this, we loaded a one-frequency shell trap at different RF-amplitudes and measured atomic transitions to untrapped states. Atomic transitions to untrapped states in this case occur at $\Omega_{\text{RF}}$ and $\omega_\text{RF}\pm \Omega_{\text{RF}}$. From here one can calibrate the RF-field amplitude.
We also determined $\Omega_1$ and $\Omega_2$ from the frequency difference between  transitions involving different number of RF photons (see Eq.\,\ref{eq:resonanceCondition}).
Initially, we set $\Omega_{1}/2\pi\approx \Omega_{2}/2\pi \approx 370$\,kHz and $\alpha=22$\,G/cm
calculated from Eqs.~(\ref{eq:matchingcondrabi},\ref{eq:matchingcondfreq}), and subsequently fine-tuned the frequency difference to compensate for the non-linear effects and imperfections in the RF-generation.
The final frequencies were $\omega_{1}/2 \pi=2.267900$\,MHz and $\omega_{2}/2 \pi=2.258910$\,MHz.
The red line is a linear fit.
Within the accuracy of the fit, the shift in angular transition frequency is equal to the shift in $\Delta \Omega$, which implies a sensitivity to the different amplitude of the driving field of 700\,kHz/G. Note that, from Fig.~\ref{fig:noisesourcesfig} b) we would expect a slightly smaller line-shift at a gradient of $\alpha=22$\,G/cm. However, in this experiment it is not clear that the traps are matched to the degree of precision that the calculation in Fig.~\ref{fig:DeltaOmegaShift} assumes. Moreover, we have not taken into account non-linear terms, which may not be negligible.
In any case, this poses strict limits on the stability of the RF generation.
For example, a coherence time of 1\,s  would require a stability of the RF magnetic field that defines $\Delta \Omega$ of, at least, $\sim 10^{-5}$.
Drifts in the gain of amplifiers or in the quality factor of resonators make such stability very hard to achieve in the case of linear polarization.
In the case of the shell traps, though,  this condition applies to the two circularly polarized fields of orthogonal polarization.
Both are generated from the the same two pairs of Helmholtz coils and differ only slightly in frequency, with the polarization determined by the relative phase of the RF in the two coils.
A change in any single coil or amplifier will cause a change in the polarization and amplitude of both RF components equally, thus considerably reducing its impact on the transition frequency.

\section{Conclusions and Outlook
\label{sec:ConclusionsAndOutlook}}

We introduced bi-chromatic adiabatic potentials as a novel system for guided/trapped atom interferometry.
We introduced a simple model for the matching condition for the two shells and demonstrated a method of fine-tuning, which results in approximately an order of magnitude reduction in the line-width of the transition between the dressed states $|1,-1\rangle$ and $|2,1\rangle$ of $^{87}$Rb compared to the single frequency case.
The resulting coherence of time of 20\,ms has been achieved for a magnetically sensitive state, without any shielding or active stabilization of either DC or RF-fields.
First calculations indicate that our set-up should be adapted for lower RF-frequencies and higher RF-amplitudes, both within reach experimentally.
We expect that with proper care the current coherence time can be extended by one or two orders of magnitude.
The demonstrated insensitivity of the transition to DC-fields at least up to 3\,G allows one to move the shell simply by applying a homogeneous field.
Applications will include an imaging atom-interferometric surface probe for DC and AC fields.
For this, one would load a BEC into a large-diameter shell-trap, approach it to a surface, and image the spatial dependence of the atom-interferometric sequence.
Further theoretical and experimental work will be required  for the shell trap to fulfill its potential, which will also serve as a testbed for TAAP-based guided atom interferometry \cite{Pandey2019HypersonicRings, Navez2016}.
The non-linear nature of the magnetic field dependence on both DC and RF suggest the existence of a \textit{magic} set of parameters for the bi-chromatic shell trap clock transition, where the dependence on the strength of the RF-fields vanish linearly.

\vspace{2cm}

\section{Acknowledgements}

We gratefully acknowledge very fruitful discussions with Thomas Fernholz, Barry Garraway and German Sinuco. This work is supported by the project ``HELLAS-CH'' (MIS 5002735) which is
implemented under the  ``Action for Strengthening Research and Innovation
Infrastructures'', funded by the Operational Programme ``Competitiveness,
Entrepreneurship and Innovation'' (NSRF 2014-2020) and co-financed by
Greece and the European Union (European Regional Development Fund). SP acknowledge financial support from the Greek Foundation for Research and
Innovation (ELIDEK) in the framework of project, \emph{Guided Matter-Wave
Interferometry} under grant agreement number 4823 and General Secretariat
for Research and Technology (GSRT). GV received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie Grant Agreement No  750017.
The author(s) would like to acknowledge the contribution of the COST Action CA16221.

\section{Contributions}

WK conceived the experimental ideas.
HM, SP, GV carried out the experiments and performed the data analysis.
HM, GV, WK developed the theory.
All authors contributed to the result discussion and manuscript writing.


\section{Introduction}
One of the central phenomena of interest in lexical tone production is the deviation of their surface realizations from canonical templates of tone categories\cite{Xu1997, Prom-on2009, Surendran2007}. In a tone language, different tone categories differing in pitch movements can distinguish different lexical meanings of a syllable (e.g., in Mandarin, the syllable ``ma" in a high level pitch contour means ``mother", whereas the same syllable spoken in a falling pitch contour means ``to scold"). Even though each tone category is defined with a general pitch contour profile (such as level, rising, falling, etc.), they typically exhibit great variability in spontaneous speech. As an example, Figure \ref{fig:tone1} shows many different realizations of Mandarin tone 1, observed during speech production experiments in the lab.

Previous works in phonology, speech prosody, and tone recognition have investigated this variability by asking questions such as: (1) What factors contribute to the variability in tone production \cite{Xu1997}?(2) How can we model the tone contour trajectory in synthesized speech \cite{Prom-on2009}? (3) What features can we use to improve the accuracy of automatic tone recognition \cite{Surendran2007}? Each of the works was driven by a particular set of theoretical or practical motivations and offered us a slice of understanding into the problem.

In this work, we are interested in looking at the tone variability problem from a data mining perspective: we explore the structure and distribution of tone contour shapes within a large amount of data. By taking a data mining approach, we contrast our work with those works that focus on tone recognition or tone learning (either by machine or by human): we seek to extract tone patterns of empirical significance from a large data set of tones from spontaneous speech.

Working with the MCPST corpus (see Section \ref{sec:data}) of Mandarin newscast speech (about 100,000 tones), we ask two questions: (1) For each tone category, what are the (coarse) types/classes of tone contour shapes we observe in this corpus? (2) For a particular tone category, what linguistic factors caused the same tone to be realized as these different types of shapes?

Inspired by works in natural language processing (NLP), we further extend these research questions in two directions. First, we extend our investigation of tone categories into a series of $n$ consecutive tones, or tone $n$-grams. $N$-grams is a classic technique in NLP language modeling\footnote{Readers may refer to the classic NLP textbook chapter if needed: https://web.stanford.edu/~jurafsky/slp3/3.pdf.}, whereas in the current context, we study tone $n$-grams due to the importance of context in tone variability \cite{Xu1997}: a tone category maybe realized differently depending on their neighboring tones. What can we learn from data mining tone contour shapes for tone unigrams, bigrams, and trigrams?

Second, to study prosody interface in MCPST data, we use automatic methods (NLP and other) to extract linguistic features from the text, including Named Entity Recognition (NER), Coreference resolution, Part-of-speech (POS) tagging, dependency parsing, and other phonological, morphological and contextual features. In order to find out the importance of these linguistic factors in shaping tone variability, we run the following experiment: given a particular tone (or tone $n$-gram) category, how well can we predict the type of tone contour shape it will take in running speech, using these linguistic features that \textit{exclude} information about the pitch contour $f0$ values?

Previous works showed that many linguistic factors (such as focus, topic, etc.) affect tone production or prosody (see Section \ref{sec:related-work}) . In this work, we extend this to a more comprehensive set of linguistic features, motivated by the information theory account of tone production. We hypothesize that there exists an information content inequality resulting from probability distribution of events in various linguistic domains (phonological, semantic, etc). These inequalities affect speakers' speech production, resulting in gradient variants of tone contour shapes in a given tone category. We investigate the relative importance of these factors in predicting the types of contour shapes any particular tone $n$-gram will take.

The rest of the paper is organized as follows. Section \ref{sec:related-work} discusses relevant previous works. Next we describe the data used in this paper in Section \ref{sec:data}. In order to characterize the types of contour shapes a tone $n$-gram will take, we develop a method to derive clusters of tone contour shape types using network analysis (Section \ref{sec:shape}).  In Section \ref{sec:features}, we discuss feature engineering and feature extraction from various linguistic domains (syntax, morphology, semantics, information structure, etc.). Section \ref{sec:ML} reports machine learning experiments and results on predicting tone contour shape types and the analysis on feature importance. Finally, in Section \ref{sec:discussion} we discuss the implications of this work in the context of information theory and phonological theory of speech and tone production.

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.1]{tone1.png}}
 \caption[Samples of Mandarin Tone1 by one speaker in lab speech]{Samples of Mandarin Tone1 by the same speaker in lab speech. Data source: \cite{Xu1997}. The canonical contours of Mandarin Tone 1,2,3,4 are: high level, low rising, low dipping, high falling, where low and high denote the pitch starting point of the tone. }
  \label{fig:tone1}

\end{figure}

\section{Related Work} \label{sec:related-work}

There has been a long line of research on the variability of tone contour shapes as well as interfacing between other linguistic factors and prosody \cite{KLi09,buring2013}. In linguistic research of Mandarin tones, most works have focused on the effect of local tonal context (e.g., neighboring tones and pitch range, such as \cite{Gauthier2007, Xu1997}) and broader context (e.g., focus, topic, information structure, long term $f_0$ variations, such as \cite{Xu2004, Liu2006, Wang2011}). The data in these works usually consisted of a small number of tone observations obtained in speech production experiments in the lab. They have informed later works on improving the performance of supervised or unsupervised tone recognition (\cite{Levow2005, Surendran2007} etc.). Other works such as \cite{Surendran2007} and \cite{KYu11} have shown the importance of signals in speech outside of $f_0$ for tone recognition and learning.

In the PENTA  \cite{Xu1997, Xu2005}  and qTA (quantitative target approximation) models \cite{Prom-on2009}, the surface $f_0$ contour is viewed as the result of asymptotic approximation to an underlying pitch target, which can be a static target (High or Low) or a dynamic target (Rise or Fall). An important contribution of the qTA is that it provides a mathematical model to account for the process of generating of a particular realization of a tone template, defined by a pitch target (with slope and intercept parameters) and the acceleration rate. As such, the specific shape of the contour then would depend on the starting pitch, ending pitch target, and how fast the pitch moves.

A fundamental theoretical question is how should we view the underlying factors that account for the tone surface variability. Previous research exhibits two opposing theories to this question. The first approach \cite{Cooper1985, Cooper1981} postulates a direct link between communicative functions and surface acoustic forms by finding the acoustic correlates of certain communicative functions, such as focus, stress, newness, questions, etc. Such approaches have met criticisms from phonologists \cite{Ladd1996, Liberman1984}, who argue that prosodic meanings are not directly mapped onto acoustic correlates. Instead, intonational meanings should be first mapped onto phonological structures, which is in turn linked to surface acoustic forms through phonetic implementation rules. In this work, we attempt to show a new middle ground between these two theories.

\section{Data} \label{sec:data}

All the data in this work comes from the Mandarin Chinese Phonetic Segmentation and Tone (MCPST) corpus \footnote{https://catalog.ldc.upenn.edu/LDC2015S05}, developed by the Linguistic Data Consortium (LDC). It contains 7,849 Mandarin Chinese newscast speech ``utterances" and their phonetic segmentation and tone labels. Utterances are defined as the time-stamped between-pause units in the transcribed news recordings. We used the auto-correlation algorithm implemented in Praat \footnote{Boersma, Paul and Weenink, David (2019). Praat: doing phonetics by computer [Computer program]. Version 6.0.48, retrieved 17 February 2019 from http://www.praat.org/} for $f_0$ (pitch) estimation from speech audio signal. We obtained $f_0$ pitch contour data for 100,161 syllables. After pre-processing the pitch tracks (e.g., speaker-dependent normalization, $f_0$ outlier detection and removal, pitch interpolation, downsampling), we generate tone unigram, bigram and trigram $f_0$ data sets, giving rise to a total of 75 unigram (5), bigram (16), and trigram data sets for the prediction task (54) \footnote{Mandarin Chinese has four regular tone categories plus one neutral tone. Since neutral tones occurs infrequently, we did not include them in the analysis of tone bigrams and trigrams due to data sparsity in the conditional distributions. Similarly, we also excluded ngrams categories where the data points are sparse.}. The total number of tone $n$-grams in these data sets are on the order of 250k. All tone unigram, bigram, and trigram $f_0$ vectors are downsampled to length of 30, 100, and 200 samples respectively.

\section{Deriving tone $n$-gram contour shape types through network analysis} \label{sec:shape}
\subsection{Problem formulation}

We define a tone $n$-gram category as a consecutive sequence of $n$ tones $t_i$, for $i=1,...,n$, where $t_i$ $\in$ $\{0,1,2,3,4\}$, the five tone categories of Mandarin. In this paper we restrict $n$ to $\{1,2,3\}$. Given the set $S$ (represented as a network) that contains all observations of $f_0$ vectors that belong to a particular tone $n$-gram category, an algorithm $A$, defined in this section, partitions $S$ into $k$ clusters, $c_1,c_2,...,c_k$, where all tone contours within $c_i$ are highly similar to each other, and members of $c_i$ maximally distinct from $c_j$ for $i \neq j$. For a particular tone $n$-gram category, we define the centroid $f_0$ vector of $c_i$ to be its tone contour shape type $t_i$.

If we denote $C$ to be set of types \{$c_1,c_2,...,c_k$\}, our goal in this section is to describe the algorithm $A$ that learns a function $g: S \rightarrow C$. We adapt an algorithm first proposed by  \cite{Gulati16}, which has been shown to be effective in identifying clusters in time-series data such as pitch contours. It also has several advantages over baseline algorithms such as $k$-means clustering, including outlier pruning and no need to determine the number of clusters before hand.

On a high level, this method represents all tone contours in a data set as a fully connected network $G$. It then filters $G$ using heuristics based on the pairwise similarity of tone contour shapes. After the filtering step, only those nodes that have a similarity score beyond a threshold will remain connected. It then leverages network community detection algorithms to optimize the community structure, effectively deriving tone contour shape types $T$ = \{$t_1,t_2,...,t_k$\}.

\subsection{Network construction}
To construct the network as described above, we first partition all data in MCPST corpus by their tone $n$-gram category. For each category, we construct a network where each node stores the $f_0$ vector of an observation, and the edge between two nodes holds the Euclidean distance between the two $f_0$ vectors as weights. We derive an undirected, weighted and fully connected network $G$ of tone $n$-gram patterns for each $n$-gram category.

\subsection{Network filtering}
In this step, we take a fully connected network $G$ of a given tone $n$-gram category and use a principled method to remove edges from the network. Our goal is to find an appropriate threshold so that all edges whose weights (distance between two tone $n$-gram $f_0$ vectors) are greater than the threshold will be cut. In the resulting network, only those nodes representing similar enough tone contour shapes will remain connected.

Specifically, we decide the threshold value by a six-step process: (1) We search for the appropriate threshold in the set of values $\Phi \in \{1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5\}$ for bigrams and trigrams, $\Phi \in \{0.2,0.4,0.6,0.8\}$ for unigrams. These values are empirically chosen; (2) we iterate over this set of values and each time apply a threshold to the network; (3) after we applied the threshold we convert the network to a unweighted network $G'$ where only those nodes that have a distance below the threshold will remain connected; (4) we produce a randomized network $G_r$ by randomly swapping edges from $G'$ $k$ times while keeping the degree of the nodes constant, where $k$ is equal to the number of edges in $G'$. This can be seen as producing a maximally random network given the degree distribution of the current network; (5) we compute the difference in Clustering Coefficient (CC) of both $G'$ and $G_r$; (6) after repeating this for all values in $T$, we pick the threshold that has the largest difference of CCs.  \footnote{Clustering coefficient (CC) measures the extent to which the nodes in a network tend to cluster together. Intuitively, it expresses how saturated the network is --- how many of the possible connections are actually expressed. The CC for a network of $k$ nodes and $n$ edges is computed as:

\[
    CC = \frac{2n}{k(k-1)}
\]

}

\subsection{Community detection}
We use the Louvain algorithm \cite{Blondel2008} to perform community detection, in order to partition the filtered network derived from last step into communities (clusters) $C_1, C_2, ..., C_k$. We pre-tuned the hyperparameters in the network filtering step so that it will result in a small number ($n<10$) of tightly connected medium-sized communities. Figure \ref{fig:num_labels} shows a histogram of number of shape classes for unigram, bigram, and trigram data sets.

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.56]{num_labels.pdf}}
 \caption{Histogram of number of shape classes in $n$-gram data sets.}
  \label{fig:num_labels}

\end{figure}

\subsection{Outlier community filtering}
We propose an extra step of outlier community filtering before deriving our final contour profile classes. In this step, we use a heuristic threshold of $t=10$ to filter out any communities (clusters) with a size less than $t$.

\subsection{Evaluation of tone contour profile classes}

We have leveraged the intrinsic structure of the tone data to derive the tone contour shape types for each tone $n$-gram category. Figure \ref{fig:panel} shows examples of learned clusters of tone contour shape types from two $n$-gram categories of tone unigram, bigram, and trigram, respectively. Without declaring any cognitive or phonological significance of these clusters, these resulting clusters should reflect the similarities of tone contour shapes within any given tone $n$-gram category: those that are highly similar are grouped into the same cluster. This is an intrinsic property derived from the above method, and is a necessary property sufficient for carrying out the subsequent experiment on predicting the tone contour shape types from linguistic factors.

Nonetheless, we propose two different ways to evaluate the validity of these clusters. First, in the following experiments, we show that we are able to predict these learned tone contour shape types significantly better than randomly assigned clusters (Section \ref{sec:results}, Figure \ref{fig:acc-box}). Second, we train a decision tree classifier to predict the shape type of a given tone $n$-gram using its $f_0$ vector and obtained a mean accuracy of 92\% (following \cite{W16-2001}). This indicates how well these tone contour shape types can be predicted \textit{with} complete information of its pitch trajectory, which will serve as an upper bound to our next prediction task using linguistic factors \textit{without} information about pitch movements $f_0$ values.

\begin{figure*}%
\subfigure[Tone 1]{%
  \label{fig:Total_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/tone1.pdf}%
}%
\hspace*{\fill}
\subfigure[Tone 2]{
  \label{fig:Num_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/tone2.pdf}%
}%
\hspace*{\fill}
\subfigure[Tone 2-3]{
  \label{fig:Raw_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/23.pdf}%
}

\subfigure[Tone 4-2]{%
  \label{fig:Num_Raw_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/42.pdf}%
}%
\hspace*{\fill}
\subfigure[Tone 3-4-2]{
  \label{fig:Total_Raw_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/342.pdf}%
}%
\hspace*{\fill}
\subfigure[Tones 2-2-4]{
  \label{fig:Num_Total_scatter}%
  \includegraphics[width=0.31\textwidth]{shapes/224.pdf}%
}

\caption{Example contour shape type clusters found in randomly selected tone $n$-gram categories, with two examples from each of tone unigram, bigram, and trigram. Each cluster is represented by a mean pitch vector with error bars (as shown in the legend). Clusters are indexed by integers shown in legend. X-axis shows number of samples (discrete time index) for the tone contour $f_0$ vector. Y-axis shows the speaker-normalized pitch values.}\label{fig:panel}
\end{figure*}

\section{Linguistic features} \label{sec:features}

After we obtained tone contour shape types emerged from each tone $n$-gram category in the corpus, we now describe the linguistic features used to predict which type of shape the tone $n$-gram will take. All syntactic and semantic features are extracted using Stanford CoreNLP for Chinese \cite{manning-EtAl:2014:P14-5}.

\subsection{Syntactic features}
Syntax and prosody has been the subject of investigation in (psycho)linguistic studies  \cite{Bard99thedissociation}.  We extract the part-of-speech (POS) tags for all syllables in a tone $n$-gram. In addition, we also extract the dependency function \cite{DBLP:conf/emnlp/ChenM14} of all syllables in the tone $n$-gram. Therefore there are $2*N$ syntactic features where N is the number of syllables included in the tone $n$-gram data under consideration. The original tag set used in CoreNLP comes from Penn Chinese Treebank\footnote{https://catalog.ldc.upenn.edu/LDC2013T21} and is too fine grained. To avoid data sparseness, we collapsed several categories. For both POS tags and dependency edge function categories, we compute their distributions using the original tag set and we collapse any categories that appear less than 5 times in the data. For POS tags we mapped the original 33 tags onto 5 categories. For dependency functions, we collapsed all tags with a subcategory separated by a colon (e.g., ``advmod:loc'', ``advmod:rcomp'', mapped to ``advmod'' etc.).

\subsection{Semantic features}
We extract two semantic features for a tone $n$-gram data point: (1) whether the tone $n$-gram includes a named entity; (2) whether the tone $n$-gram includes a singleton (as opposed to being part of a coreference chain in the discourse). Semantic features such as information structure have been postulated to have an effect on the prosody domain \cite{buring2013}. In particular, given information may encode prosodic features different from new information. This could also apply to named entities vs. non-named entities. Named entities points to definite, specific objects in the real world. Whether the token is a singleton (i.e., does not co-refer to an entity with another mention in the text) or part of a coreference chain can be correlated with information structure \cite{recasens-etal-2013-life}. That is, a singleton may signify new information in discourse, while a non-singleton is part of a coreference chain with potential antecedent or anaphor, pointing to potentially a different information structure. Both can have distinguishing effects on the mental representations and the production of speech prosody (indirectly related to redundancy in \cite{Aylett2004}).

\begin{table*}
 \caption[Feature set overview]{Feature set overview. \textmd{1...N indicates this feature is computed for all syllables in the $N$-gram. Total number is N*10+7 features, 37 for trigrams and 27 for bigrams, etc.}}
\small

\begin{center}
\begin{tabular}{ |c|c|c|c|c| }

 \hline
Syntactic & Morphological & Semantic & Phonological & Others\\
 \hline
\hline
$POS\_Tag_{1...N}$ & $Tok\_Bound_{1...N}$ & is\_entity & $is\_nasal_{1...N}$ & sent\_position \\
$Dep\_Func_{1...N}$ & & is\_singleton & $is\_dipthong_{1...N}$ & start\_pitch \\
 & & & $is\_round_{1...N}$ & end\_pitch\\
 & & & $is\_front_{1...N}$ & prev\_tone\\
  & & & $is\_back_{1...N}$ & next\_tone\\
   & & & $is\_high_{1...N}$ & \\
    & & & $is\_low_{1...N}$ & \\

\hline
\end{tabular}

  \label{tab:trigramFeatures}
\end{center}
\end{table*}

\subsection{Morphological features}
In Mandarin Chinese, each word usually consists of one to four syllables. Building on the intuition that the first syllable is usually spoken with higher prominence (e.g., neutral tone, which does not carry stress, only occurs on word-final positions), we extract morphological features for each syllable in the given tone $n$-gram: whether they cross word boundary or not. There are $n$ features in this category in total.

\subsection{Phonological features}
A basic representation of phonological features is the identity of phonemes in each syllable of the $n$-gram. However, due to the sparseness of this feature representation, we have designed 7 binary features to encode the phonological properties of the syllables in the tone $n$-gram: (1) whether the syllable includes a nasal; (2)  whether the syllable includes a dipthong; (3) whether the syllable includes a high vowel; (4)  whether the syllable includes a low vowel; (5)  whether the syllable includes a front vowel; (6)  whether the syllable includes a back vowel; (7)  whether the syllable includes a round vowel. In addition, we add two contextual tone features: the tone identity of the previous and following syllables of the tone $n$-gram in question.

\subsection{Other features}
We add two pitch features to the feature set: the beginning and ending pitch of the tone $n$-gram. This is based on the notion in generative Mandarin tone modeling (Parallel ENcoding and Target Approximation model, or PENTA) that in speech production, the actual realized tone shape of a given tone category highly depends on the starting point of the pitch contour and its distance to the actual pitch target of the current tone, which affects its course of trajectory when it approximates the target \cite{Prom-on2009}. An additional feature to be included is the position of the current tone $n$-gram within the context of the current sentence as a percentage. It is a known effect that pitch tends to downdrift in speech production as sentence progresses \cite{Wang2011}. Therefore, we also want to account for the effect of sentence position.

\subsection{Bag of features}
In this task, we note that the unit of feature extraction is not as straightforward as it would be in classic NLP tasks. That is, instead of a typical syntactic constituent (word, phrase, sentences) as the feature extraction unit, here, our target is tone $n$-gram, a sequence of $n$ syllables that may or may not be a syntactic constituent. As described above, in many features we have adopted a  ``Bag of features'' approach (similar to the speech coreference resolution work in Roesiger et al.  \cite{Roesiger2015}) where each feature describes whether the $n$-gram contains a certain target value in any position. For some other features, we simply use a set of $n$ features applied to each syllable in the $n$-gram in question. These are precisely described in Table \ref{tab:trigramFeatures}.

\section{Predicting tone contour profiles} \label{sec:ML}

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.58]{new_results.pdf}}
 \caption{Distribution of classification accuracies for all 75 data sets in unigram(U), bigram(B), trigram(T). For each $n$-gram the labels on the x-axis are: data:full feature set; dfp:start/end pitch only baseline; no\_syn:without syntactic features; no\_Ntone:without prev/next tone features; no\_pitch:without start/end pitch features; random:baseline with randomly assigned labels; MLE:MLE baseline.}
  \label{fig:acc-box}

\end{figure}

\subsection{Experimental setup}

For each of the 5 unigram, 16 bigram, and 54 trigram data sets, the extracted linguistic feature vector ($f_1,f_2,...,f_m$) forms the input space $X$. The contour shape types $T$ forms the output space $Y$. Our goal is to learn a function $h: X \rightarrow Y$ minimizing expected loss. We use SVM with linear kernel so that we can extract feature importance in subsequent analyses. Each data set is randomly split into 90\slash10 for train and test. Since the classes are balanced in $Y$, we evaluate the classifier performance directly with accuracy on the test sets.

\subsection{Results} \label{sec:results}
To visualize model performances on a large number of tone $n$-gram data sets (75), we choose the boxplot because of its efficiency to convey statistical information of the distribution of the results across all data sets. Figure \ref{fig:acc-box} gives an overview of our proposed model classification accuracies for all unigram(U), bigram(B), and trigram(T) data sets, as compared to several baselines. In this figure, \texttt{data}, \texttt{dfp}, \texttt{no\_syn}, \texttt{no\_Ntone}, and \texttt{no\_pitch} denote results using different sets of features: full set, start/end pitch only, no syntactic features, no prev/next tone features, and no start/end pitch features, respectively. The \texttt{random} baseline uses the same set of linguistic features as our proposed model but the target tone contour shape type is randomly assigned (while keeping the number of shape types constant). Finally, the Maximum Likelihood Estimation (\texttt{MLE}) baseline reflects chance level performance if linguistic factors are independent from the output tone contour shape types, and is calculated as $1/d$, where $d$ is the number of output classes in a data set.

First, the proposed model significantly outperforms MLE baselines for unigrams, bigrams, and trigrams. Second, the accuracy for predicting learned tone contour shape types is significantly higher than randomly assigned clusters. This serves as a sanity check for the validity of these learned tone contour shape types. Overall, this result supports the hypothesis that a variety of linguistic and contextual features are strongly correlated with the realization of a particular category of tone $n$-grams. In particular, we observe that using the set of features excluding the syntactic features (POS tags and dependency functions) allows the model to achieve the best median accuracy in all $n$-grams (\texttt{no\_syn} baseline). In Section \ref{sec:feature-imp}, we provide a more detailed analysis and discussion of feature importance.

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.56]{unigram_coef_new.pdf}}
 \caption[Feature weights for all unigram data sets]{Feature weights for all unigram data sets.}
  \label{fig:unigramcoef}

\end{figure}

A significant trend to note in Figure \ref{fig:acc-box} is that we observe a negative correlation between the model accuracy and the baseline (MLE, random) accuracy as the value of $n$ becomes larger in $n$-grams. This is striking because it indicates a decrease in predictive power as $n$ grows larger. Moreover, the variance on model accuracy also increases as $n$ becomes larger. From inspecting the tone contour shape types we obtained (such as those showed in Figure \ref{fig:panel}), we attribute this to three factors from the data perspective: (1) the dimensionality of unigram, bigram, and trigram $f_0$ vectors in our data set are different (increase); (2) the number of data sets also increases as $n$ grew larger; (3) the complexity of tone contour shapes tend to increase from unigram to bigram to trigram.  On a linguistic level, we hypothesize that the longer the window of $n$-grams, the stronger an effect of unaccounted factors come into play (longer range prosodic factors such as focus and topic, as demonstrated in \cite{Xu2004}).

\subsection{Feature importance} \label{sec:feature-imp}

To analyze feature importance, we extracted weight vectors (coefficients) associated with all features in the linear SVM classifiers (following \cite{Guyon2002}) for all data sets in each of the $n$-grams. We take the absolute values of feature weights and normalize them to be comparable across data sets. We then aggregated feature weights for all classes and across all data sets in a given $n$-gram.

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.56]{bigram_coef_new.pdf}}
 \caption[Feature weights for all bigram data sets]{Feature weights for all bigram data sets.}
  \label{fig:bigramcoef}

\end{figure}
Figure \ref{fig:unigramcoef}, \ref{fig:bigramcoef}, \ref{fig:trigramcoef} show the distribution of feature weights (importance) for all features for unigram, brigram, and trigram data. Since feature importance values tend to behave similarly according to their linguistic domains, we group them together instead of analyzing the features individually. We observe three levels of feature importance based on the weights, consistent among different $n$-grams: (1) High: starting pitch, ending pitch, and sentence position (especially when $n>1$). The importance of starting and ending pitch is consistent with the qTA model of Mandarin tones \cite{Prom-on2009}. The latter (sentence position) is consistent with the effect of downdrift \cite{Wang2011}. (2) Medium: Phonological features and morphological token boundary features, as well as coreference (singleton) and entity, which are information structure of the discourse and semantics (givenness, newness of information in speech). (3) Low: Syntactic (pos tag, dependency functions) and contextual features (previous and next tone). This is consistent with \cite{Bard99thedissociation} and \cite{Surendran2007}\footnote{Specifically, \cite{Bard99thedissociation} showed the dissociation of syntax with de-accent in spontaneous speech, and  \cite{Surendran2007} showed that context did not help tone recognition as much as expected.}.

To have a more detailed understanding of the feature importance, Figure \ref{fig:acc-box} shows how the model performs with partially ablated feature sets across different $n$-grams. First, to understand the role of start/end pitch vs. non-pitch linguistic features, we observe that the \texttt{dfp} baseline (using only start/end pitch features) has lower results than other models with linguistic features. This is true for all the markers on the boxplots (min, max, first quantile, median, third quantile) when comparing \texttt{dfp} to \texttt{no\_syn}. However, there is also considerable overlap in these accuracy distributions to different degrees as $n$ varies, which indicates cases where the \texttt{dfp} baseline outperforms the other models with more linguistic features. To see this possibility, we plotted the difference (delta) in accuracies values for all data sets for the pair of baselines \texttt{no\_syn} - \texttt{dfp} in Figure \ref{fig:diffs}. It shows that for 80\% of the data sets, the \texttt{no\_syn} baseline with linguistic features outperforms the pitch-only baseline in most data sets by a margin of less than 20\% in accuracy improvements.

Second, comparing across different $n$-grams, the \texttt{no\_pitch} baseline (only linguistic features) performs worse in unigram, and the best in trigrams. This shows that the pitch feature is less important when $n$ becomes larger in $n$-grams. The same trend is observed in the weights of the pitch features in feature importance. This observation is also consistent with \cite{Prom-on2009}: the target approximation of tones only operates on the syllable units. Therefore the effect of start/end pitch should diminish when $n\textgreater1$. Third, even though feature weights are small for previous/next tones, in these results we didn't see an improvement when we exclude these features in the \texttt{no\_Ntone} baseline. Finally, we confirm the importance of the non-pitch linguistic features since the \texttt{no\_pitch} baseline significantly outperforms the \texttt{random} and the \texttt{MLE} baselines.

\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.56]{trigram_coef_new.pdf}}
 \caption[Feature weights for all trigram data sets]{Feature weights for all trigram data sets.}
  \label{fig:trigramcoef}

\end{figure}

\section{Discussion} \label{sec:discussion}

In this work, we first described a method to mine tone contour shape types from a large amount of tone data. These shape types, as well as those mined from tone $n$-grams with larger values of $n$ ($n\textgreater3$), can be further analyzed with linguistic knowledge to better understand the behaviors of tones in different tonal contexts in future works. We also showed that using linguistic and contextual factors, we can predict with reasonable accuracy the contour shape type that a tone $n$-gram will take in the MCPST data. We analyzed the feature importance to form a relative ranking of the linguistic factors. These results should be interpreted with caution because they are constrained by the particular representations of the linguistic features used in this study, as well as the accuracy of the NLP softwares used to extract them. Nonetheless, we envision that by mining correlates between speech prosody and automatic analysis of linguistic features extracted from the data, this line of work could have potential applications in improving the quality and naturalness of prosody in speech synthesis such as Text-To-Speech (TTS) technologies.

Previous works targeting information theory and information structure in prosody domain have largely looked at acoustic correlates directly, such as accent and duration, all of which may in turn have an impact on the shape of tone contours in speech production. Therefore, looking at tone contour shapes can be thought of as a different level of manifestation of such phenomena, an amalgamation of single dimension acoustic correlates (e.g., duration and intensity). It is also a level that is more difficult to quantify and measure in the traditional linguistic/phonetic investigations on a smaller scale. One possible extension of this work in the future is to look at specific ways tone contour shapes correlate with particular features when they are perturbed in a certain direction. Moreover, it is also of interest to demonstrate this change in a quantified manner using the information theory formulation.

In Section \ref{sec:related-work} we raised a fundamental theoretical question of whether there is a direct link between communicative functions and surface acoustic forms, a question that we found disagreement in literature (as summarized in \cite{Xu2005} and in Section \ref{sec:related-work} of the current paper). In this paper, we showed that by taking a data driven approach, we can predict the contour shape type of a prosodic category (such as a tone $n$-gram) using linguistic factors, even though we are not predicting its exact shape. In doing so, we give an approximate solution for a middle ground between the two theories.
\begin{figure}

\small

 \centerline{
 \includegraphics[scale=0.56]{diffs.pdf}}
 \caption{Delta (differences in accuracies): no\_syn - dfp baselines. A point above the y=0 line indicates the model with linguistic features outperforms the model with pitch only features. }
  \label{fig:diffs}

\end{figure}
\section*{Acknowledgements}
We thank Sankalp Gulati, Xavier Serra, Amir Zeldes, Elizabeth Zsiga, George Wilson, Richard Wright and Gina-Anne Levow for insights and discussions, as well as the anonymous reviewers for valuable comments on the previous versions of this paper.