We present an inference of the nuclear symmetry energy magnitude $J$, the slope $L$, and the curvature ${K}_{\mathrm{sym}}$ from combining neutron skin data on calcium, lead and tin isotopes, and our best theoretical information about pure neutron matter. A Bayesian framework is used to consistently incorporate prior knowledge of the pure neutron matter equation of state from chiral effective field theory calculations. Neutron skins are modeled in a fully quantum Skyrme-Hartree-Fock approach using an extended Skyrme energy-density functional which allows for independent variation of $J$, $L$, and ${K}_{\mathrm{sym}}$ without affecting the symmetric nuclear matter equation of state. The effect of using neutron skin data obtained with different physical probes is quantified. We argue that, given the existing data, combining the errors in quadrature is the more appropriate way to obtain unified errors for each nuclide, and in doing so we obtain 95% credible values of $J=31.3\begin{array}{c}+4.2\\ \ensuremath{-}5.9\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$, $L=40\begin{array}{c}+34\\ \ensuremath{-}26\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$, and ${K}_{\ensuremath{\tau}}=L\ensuremath{-}6{K}_{\mathrm{sym}}=\ensuremath{-}444\begin{array}{c}+100\\ \ensuremath{-}84\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ using uninformative priors in $J$, $L$, and ${K}_{\mathrm{sym}}$, and $J=31.9\begin{array}{c}+1.3\\ \ensuremath{-}1.3\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$, $L=37\begin{array}{c}+9\\ \ensuremath{-}8\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$, and ${K}_{\ensuremath{\tau}}=\ensuremath{-}480\begin{array}{c}+25\\ \ensuremath{-}26\end{array}\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ using pure neutron matter (PNM) priors. We also show that the nonpositive correlation between $J$ and $L$ induced by neutron skin data is consistent with the nuclear droplet model. Neutron skin data alone are shown to place limits on the symmetry energy parameters as stringent as those obtained from chiral effective field theory alone, and when combined the 95% credible intervals are reduced by a factor of 4--5. It is also shown that the majority of nuclear interactions used in the literature have subsaturation density dependencies that are inconsistent with the combination of PNM priors and neutron skin data. We show measurements of lead and calcium neutron skins from upcoming parity-violating electron scattering experiments at Jefferson Lab and Mainz Superconducting Accelerator should obtain total error ranges $\mathrm{\ensuremath{\Delta}}L\ensuremath{\approx}50\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ and $\mathrm{\ensuremath{\Delta}}{K}_{\ensuremath{\tau}}\ensuremath{\approx}240\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ for uninformative priors and $\mathrm{\ensuremath{\Delta}}L\ensuremath{\approx}30\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ and $\mathrm{\ensuremath{\Delta}}{K}_{\ensuremath{\tau}}\ensuremath{\approx}100\phantom{\rule{0.28em}{0ex}}\mathrm{MeV}$ for PNM priors at 67% credible bounds. Ahead of those experiments, we make predictions based on existing data on neutron skins of tin alone for the neutron skins of calcium and lead of $0.166\ifmmode\pm\else\textpm\fi{}0.008$ fm and $0.169\ifmmode\pm\else\textpm\fi{}0.014$ fm, respectively, using uninformative priors and $0.167\ifmmode\pm\else\textpm\fi{}0.008$ fm and $0.172\ifmmode\pm\else\textpm\fi{}0.015$ fm, respectively, using PNM priors.