  PICO is a concept for a NASA probe-scale mission aiming to detect or
constrain the tensor to scalar ratio $r$, a parameter that quantifies the
amplitude of inflationary gravity waves. We carry out map-based component
separation on simulations with five foreground models and input $r$ values
$r_{in}=0$ and $r_{in} = 0.003$. We forecast $r$ determinations using a
Gaussian likelihood assuming either no delensing or a residual lensing factor
$A_{\rm lens}$ = 27%. By implementing the first full-sky, post
component-separation, map-domain delensing, we show that PICO should be able to
achieve $A_{\rm lens}$ = 22% - 24%. For four of the five foreground models we
find that PICO would be able to set the constraints $r < 1.3 \times 10^{-4}
\,\, \mbox{to} \,\, r <2.7 \times 10^{-4}\, (95\%)$ if $r_{in}=0$, the
strongest constraints of any foreseeable instrument. For these models,
$r=0.003$ is recovered with confidence levels between $18\sigma$ and
$27\sigma$. We find weaker, and in some cases significantly biased, upper
limits when removing few low or high frequency bands. The fifth model gives a
$3\sigma$ detection when $r_{in}=0$ and a $3\sigma$ bias with $r_{in} = 0.003$.
However, by correlating $r$ determinations from many small 2.5% sky areas with
the mission's 555 GHz data we identify and mitigate the bias. This analysis
underscores the importance of large sky coverage. We show that when only low
multipoles $\ell \leq 12$ are used, the non-Gaussian shape of the true
likelihood gives uncertainties that are on average 30% larger than a Gaussian
approximation.
