By using monthly mean water levels at 55 sites around the Great Lakes, a regional model of vertical crustal motion was computed for the region. In comparison with previous similar studies over the Great Lakes, 15 additional gauge sites, data from all seasons instead of the 4 summer months, and 8 additional years of data were used. All monthly water levels available between 1860 and 2000, as published by the U.S. National Ocean Survey and the Canadian Hydrographic Service, were used. For each lake basin, the vertical velocities of the gauge sites relative to each other were simultaneously computed, using the least-squares adjustment technique. Our algorithm solves for and removes a monthly bias common to all sites, as well as site-specific biases. It also properly weighs the input water levels, resulting in a realistic estimation of the uncertainties in tilting parameters. The relative velocities obtained for each lake were then combined to obtain relative velocities over the entire Great Lakes region. Finally, the gradient of the relative rates for the regional model was found to agree best with the ICE-3G global isostatic model of Tushingham and Peltier, whereas the ICE-4G gradients were too small around the Great Lakes.