  We consider a harmonically trapped dilute $N$-boson system described by a
low-energy Hamiltonian with pairwise interactions. We determine the condensate
fraction, defined in terms of the largest occupation number, of the
weakly-interacting $N$-boson system ($N \ge 2$) by employing a perturbative
treatment within the framework of second quantization. The one-body density
matrix and the corresponding occupation numbers are compared with those
obtained by solving the two-body problem with zero-range interactions exactly.
Our expressions are also compared with high precision {\em{ab initio}}
calculations for Bose gases with $N=2-4$ that interact through finite-range
two-body model potentials. Non-universal corrections are identified to enter at
subleading order, confirming that different low-energy Hamiltonians,
constructed to yield the same energy, may yield different occupation numbers.
Lastly, we consider the strongly-interacting three-boson system under
spherically symmetric harmonic confinement and determine its occupation numbers
as a function of the three-body "Efimov parameter".
