The ``Brownian bees" model describes an ensemble of $N$ independent branching Brownian particles. When a particle branches into two particles, the particle farthest from the origin is eliminated so as to keep a constant number of particles. In the limit of $N\to \infty$, the spatial density of the particles is governed by the solution of a free boundary problem for a reaction-diffusion equation. At long times the particle density approaches a spherically symmetric steady state solution with a compact support. Here we study fluctuations of the ``swarm of bees" due to the random character of the branching Brownian motion in the limit of large but finite $N$. We consider a one-dimensional setting and focus on two fluctuating quantities: the swarm center of mass $X(t)$ and the swarm radius $\ell(t)$. Linearizing a pertinent Langevin equation around the deterministic steady state solution, we calculate the two-time covariances of $X(t)$ and $\ell(t)$. The variance of $X(t)$ directly follows from the covariance of $X(t)$, and it scales as $1/N$ as to be expected from the law of large numbers. The variance of $\ell(t)$ behaves differently: it exhibits an anomalous scaling $\ln N/N$. This anomaly appears because all spatial scales, including a narrow region near the edges of the swarm where only a few particles are present, give a significant contribution to the variance. We argue that the variance of $\ell(t)$ can be obtained from the covariance of $\ell(t)$ by introducing a cutoff at the microscopic time $1/N$ where the continuum Langevin description breaks down. Our theoretical predictions are in good agreement with Monte-Carlo simulations of the microscopic model. Generalizations to higher dimensions are briefly discussed.