To simulate, predict and control Electric Power Systems (EPS), engineers need tools to model the components of this highly complex interconnected network. Many efforts over the past century were dedicated to develop mathematical models for generators, transmission lines, reactive power compensators, transformers and so on. The main components of the power systems are precisely represented by mathematical models, but the loads are still a source of uncertainty in the simulations, due to their random characteristics. It is well known that conservative load models super estimate power response to voltage deviations, and, on the other hand, over-optimistic load models may underestimate stability margins, leading a system to operate too close to its limit. It is necessary to stablish load representations as close to reality as possible, in order to fully exploit grid resources. This work provides a methodology for load modeling, investigating and summarizing the steps of the process, whose can be implemented in a wide variety of ways. Data treatment, the choice of a load model representation and their parameters estimation are presented through real case studies, both for dynamic simulation and a steady state application. It is discussed how optimization and statistical inference concepts can be effective tools to reach better approximations on how load will respond to disturbances caused by voltage variations, whether these were spontaneous, due to control actions, or caused by short-circuits.