  Recently, the use of neural quantum states for describing the ground state of
many- and few-body problems has been gaining popularity because of their high
expressivity and ability to handle intractably large Hilbert spaces. In
particular, methods based on variational Monte Carlo have proven to be
successful in describing the physics of bosonic systems such as the
Bose-Hubbard (BH) model. However, this technique has not been systematically
tested on the parameter space of the BH model, particularly at the boundary
between the Mott insulator and superfluid phases. In this work, we evaluate the
capabilities of variational Monte Carlo with a trial wavefunction given by a
Restricted Boltzmann Machine to reproduce the quantum ground state of the BH
model on several points of its parameter space. To benchmark the technique, we
compare its results to the ground state found through exact diagonalization for
small one-dimensional chains. In general, we find that the learned ground state
correctly estimates many observables, reproducing to a high degree the phase
diagram for the first Mott lobe and part of the second one. However, we find
that the technique is challenged whenever the system transitions between
excitation manifolds, as the ground state is not learned correctly at these
boundaries. We improve the quality of the results produced by the technique by
proposing a method to discard noisy probabilities learned in the ground state.
