The method of simulated scores (MSS) is presented for estimating limited dependent variables models (LDV) with flexible correlation structure in the unobservables. We propose simulators that are continuous in the unknown parameter vectors, and hence standard optimization methods can be used to compute the MSS estimators that employ these simulators. The first continuous method relies on a recursive conditioning of the multivariate normal density through a Cholesky triangularization of its variance-covariance matrix. The second method combines results about the conditions of the multivariate normal distribution with Gibbs resampling techniques. We establish consistency and asymptotic normality of the MSS estimators and derive suitable rates at which the number of simulations must rise if biased simulators are used.