Disclaimer: I am assuming that you are implementing a classical Monte Carlo estimator.
The Problem
Discarding samples will change PDF of your sampling technique. You are cutting off part of the sampled domain where PDF is non-zero, which effectively leads to a trimmed version of the original PDF but implicitly re-normalized so the remaining part integrates to 1. If you don't adjust the directly evaluated PDF accordingly, it will lead to a biased estimator.
Practically speaking, implicit re-normalization increases the actual sampling PDF $p^{new}$, and if you use $p^{old}<p^{new}$ in your computations instead, the resulting Monte Carlo estimator will yield brighter values than it should:
$$ \frac{f(x)}{p^{old}(x)} > \frac{f(x)}{p^{new}(x)} $$
A solution
Since adjusting the PDF isn't an easy thing to do, you will most likely need to treat the under-surface samples as valid but with zero contribution. Whether you handle this within your BRDF/BSDF or elsewhere in the renderer is your design decision.
Zero-contribution samples will, obviously, introduce some inefficiency in your renderer.
A better solution
You can improve the efficiency of your estimator by using a better sampling technique which tries to avoid creating samples under the surface. In case of the GGX normal distribution, there have been proposed some solutions by Eric Heitz | Eugene d’Eon in paper Importance Sampling Microfacet-Based BSDFs using the Distribution of Visible Normals. I believe there was one improved version of this technique (other that the one mentioned in the "Related Work" section) but I cannot recall the name of it. Maybe someone else can add it here...