Context. We assess statistical inversion of asteroid rotation periods, pole orientations, shapes, and phase curve parameters from photometric lightcurve observations, here sparse data from the ESA Gaia space mission (Data Release 2) or dense and sparse data from ground-based observing programs. Aims. Assuming general convex shapes, we develop inverse methods for characterizing the Bayesian a posteriori probability density of the parameters (unknowns). We consider both random and systematic uncertainties (errors) in the observations, and assign weights to the observations with the help of Bayesian a priori probability densities. Methods. For general convex shapes comprising large numbers of parameters, we developed a Markov-chain Monte Carlo sampler (MCMC) with a novel proposal probability density function based on the simulation of virtual observations giving rise to virtual least-squares solutions. We utilized these least-squares solutions to construct a proposal probability density for MCMC sampling. For inverse methods involving triaxial ellipsoids, we update the uncertainty model for the observations. Results. We demonstrate the utilization of the inverse methods for three asteroids with Gaia photometry from Data Release 2: (21) Lutetia, (26) Proserpina, and (585) Bilkis. First, we validated the convex inverse methods using the combined ground-based and Gaia data for Lutetia, arriving at rotation and shape models in agreement with those derived with the help of Rosetta space mission data. Second, we applied the convex inverse methods to Proserpina and Bilkis, illustrating the potential of the Gaia photometry for setting constraints on asteroid light scattering as a function of the phase angle (the Sun-object-observer angle). Third, with the help of triaxial ellipsoid inversion as applied to Gaia photometry only, we provide additional proof that the absolute Gaia photometry alone can yield meaningful photometric slope parameters. Fourth, for (585) Bilkis, we report, with 1-σ uncertainties, a refined rotation period of (8.5750559 ± 0.0000026) h, pole longitude of 320.6° ± 1.2°, pole latitude of - 25.6° ± 1.7°, and the first shape model and its uncertainties from convex inversion. Conclusions. We conclude that the inverse methods provide realistic uncertainty estimators for the lightcurve inversion problem and that the Gaia photometry can provide an asteroid taxonomy based on the phase curves.