This paper studies numerical methods for the approximation of elliptic PDEs with lognormal coefficients of the form
−
div
(
a
∇
u
)
=
f
-\operatorname {div}(a\nabla u)=f
where
a
=
exp
(
b
)
a=\exp (b)
and
b
b
is a Gaussian random field. The approximant of the solution
u
u
is an
n
n
-term polynomial expansion in the scalar Gaussian random variables that parametrize
b
b
. We present a general convergence analysis of weighted least-squares approximants for smooth and arbitrarily rough random field, using a suitable random design, for which we prove optimality in the following sense: their convergence rate matches exactly or closely the rate that has been established in Bachmayr, Cohen, DeVore, and Migliorati [ESAIM Math. Model. Numer. Anal. 51 (2017), pp. 341–363] for best
n
n
-term approximation by Hermite polynomials, under the same minimial assumptions on the Gaussian random field. This is in contrast with the current state of the art results for the stochastic Galerkin method that suffers the lack of coercivity due to the lognormal nature of the diffusion field. Numerical tests with
b
b
as the Brownian bridge confirm our theoretical findings.