The oligomerization of lysozyme in aqueous solution was investigated by Monte Carlo simulations as a function of protein concentration, pH, and electrolyte screening. Lysozyme was modeled as a hard sphere with embedded pH-dependent discrete charges and with an attractive 1/r(6)-potential representing nonspecific short-range attraction. The magnitude of the 1/r(6)-potential was adjusted to reproduce experimental second virial coefficients. Radial distribution functions, structure factors, cluster size distributions, and orientation correlations were determined at various conditions. It was observed that increasing protein concentration, or decreasing the electrostatic repulsion between protein molecules by either reducing the protein charge or increasing the ionic strength, promoted cluster formation. Structure factors and equilibrium constants obtained were compared to those obtained experimentally and were found to capture the experimentally obtained effects of pH and ionic strength. The influence of the location of the hydrophobic site was also examined.