= Solution
Use the <periodic covariance function>
$$
k_P(t,t')=A^2\exp\left[-\frac2{\ell^2}\sin^2\left(\frac{\pi(t-t')}P\right)\right],\qquad\mathcal H=(A,\ell,P),\quad A,\ell,P>0.
$$
Its positive semidefiniteness follows by pulling back a <Gaussian kernel> to the circle $(\cos(2\pi t/P),\sin(2\pi t/P))$. A zero <Gaussian process> mean is an ensemble statement, whereas an exactly zero long-term average for every periodic <light curve> is stronger. To impose the latter literally, use a <periodic Gaussian process with zero period average>, replacing the kernel by
$$
k_0(t,t')=k_P(t,t')-c_0,\qquad c_0=A^2e^{-\ell^{-2}}I_0(\ell^{-2}).
$$
Here $I_0$ is a <modified Bessel function>, and $c_0$ is the period average of the kernel. Subtracting each process's random period average produces this <covariance>, so it remains positive semidefinite. We use $k_0$ to enforce the stated mean subtraction; using $k_P$ gives the usual ensemble-zero-mean formulation and the same likelihood derivation.
Let $K_{ij}=k_0(t_i,t_j)$ and $R=\operatorname{diag}(\sigma_1^2,\ldots,\sigma_N^2)$. The prior latent vector is $f\sim N(0,K)$ and the independent <Gaussian noise> is $\varepsilon\sim N(0,R)$, so $y=f+\varepsilon\sim N(0,V)$ with $V=K+R$. The <Gaussian-process marginal likelihood> is
$$
\boxed{p(y\mid t,\mathcal H)=(2\pi)^{-N/2}|V|^{-1/2}\exp(-\tfrac12y^TV^{-1}y)}.
$$
Positive measurement <variances> make $V$ positive definite even if $K$ is singular. Maximize $\ell(\mathcal H)=-\frac12y^TV^{-1}y-\frac12\log|V|-\frac N2\log(2\pi)$ over all three <hyperparameters>, using a broad search in $P$ to find competing aliases, followed by joint optimization. <Cholesky decomposition> evaluates the likelihood without explicitly forming an inverse.
For a well-resolved interior maximum, let $J=-\nabla^2\ell(\widehat{\mathcal H})$ be the <observed Fisher information> in the coordinates $(A,\ell,P)$. Then $\boxed{\widehat P=(\widehat{\mathcal H})_P,\qquad\sigma_P\approx\sqrt{(J^{-1})_{PP}}}$. Equivalently, the profile-likelihood curvature in $P$ gives the same local nuisance-adjusted <variance>; keeping the other <hyperparameters> fixed would generally underestimate it. A profile interval with $\ell_{\rm prof}(\widehat P)-\ell_{\rm prof}(P)=1/2$ is a local $1\sigma$ approximation. If the period likelihood is multimodal or poorly constrained, fit the complete posterior with proper <hyperparameter> priors and report its modes and marginal credible intervals rather than a single Hessian error bar. The sampling pattern can create period aliases that a single local search misses.
Back to article page
