  Consider an $N\times N$ Toeplitz matrix $T_N$ with symbol ${a }(\lambda) :=
\sum_{\ell=-d_2}^{d_1} a_\ell \lambda^\ell$, perturbed by an additive noise
matrix $N^{-\gamma} E_N$, where the entries of $E_N$ are centered i.i.d.~random
variables of unit variance and $\gamma>1/2$. It is known that the empirical
measure of eigenvalues of the perturbed matrix converges weakly, as
$N\to\infty$, to the law of ${a}(U)$, where $U$ is distributed uniformly on
$\mathbb{S}^1$. In this paper, we consider the outliers, i.e. eigenvalues that
are at a positive ($N$-independent) distance from ${a}(\mathbb{S}^1)$. We prove
that there are no outliers outside ${\rm spec} \, T({a})$, the spectrum of the
limiting Toeplitz operator, with probability approaching one, as $N \to
\infty$. {In contrast,} in ${\rm spec}\, T({a})\setminus {a}({\mathbb S}^1)$
the process of outliers converges to the point process described by the zero
set of certain random {analytic} functions. The limiting random {analytic}
functions can be expressed as linear combinations of the determinants of finite
sub-matrices of an infinite dimensional matrix, whose entries are i.i.d.~having
the same law as that of $E_N$. The coefficients in the linear combination
depend on the roots of the polynomial $P_{z, {a}}(\lambda):= ({a}(\lambda)
-z)\lambda^{d_2}=0$ and semi-standard Young Tableaux with shapes determined by
the number of roots of $P_{z,{a}}(\lambda)=0$ that are greater than one in
moduli.
