  $\mathit{Herschel}$ extragalactic surveys offer a unique opportunity to
efficiently select a significant number of rare and massive dusty objects, and
thus gain insight into the prodigious star-forming activity that takes place in
the very distant Universe. To search for $z\geq4$ dusty star-forming galaxies,
in this work we consider red SPIRE objects with fluxes rising from 250 $\mu$m
to $500\:\mu$m (so-called "500 $\mu$m-risers"). We aim to implement a novel
method to obtain a statistical sample of "500 $\mu$m-risers" and fully evaluate
our selection inspecting different models of galaxy evolution. We consider one
of the largest and deepest ${\it Herschel}$ surveys, the Herschel Virgo Cluster
Survey. We develop a novel selection algorithm which links the source
extraction and spectral energy distribution fitting. We select 133 "500
$\mu$m-risers" over 55 deg$^{2}$, imposing the criteria:
$S_{500}>S_{350}>S_{250}$, $S_{250}>13.2$ mJy and $S_{500}>$30 mJy.
Differential number counts are in a fairly good agreement with models,
displaying better match than other existing samples. In order to interpret the
statistical properties of selected sources, which has been proven as a very
challenging task due the complexity of observed artefacts, we make end-to-end
simulations including physical clustering and lensing. The estimated fraction
of strongly lensed sources is $24^{+6}_{-5}\%$ based on models. We present the
faintest known statistical sample of "500 $\mu$m-risers" and show that noise
and strong lensing have crucial impact on measured counts and redshift
distribution of selected sources. We estimate the flux-corrected star formation
rate density at $4<z<5$ with the "500 $\mu$m-risers" and found it close to the
total value measured in far-infrared. It indicates that colour selection is not
a limiting effect to search for the most massive, dusty $z>4$ sources.
