The standard Gaussian distribution for incremental frequency data requires a constant variance which is independent of the mean. We develop a more general and appropriate method based on the Poisson distribution, which assumes different unknown variances for the frequencies, equal to the means. We explicitly include empty bins, and our method is quite insensitive to the choice of bin width. We develop a maximum likelihood technique that minimizes bias in the curve fits, and penalizes additional free parameters by objective information criteria. Various data sets are used to test three different physical models that have been suggested for the density distribution: the power law; the double power law; and the gamma distribution. For the CMT catalog of global earthquakes, two peaks in the posterior distribution are observed at moment magnitudes m*=6.4 and 6.9 implying a bimodal distribution of seismogenic depth at around 15 and 30 km, respectively. A similar break at a characteristic length of 60 km or so is observed in moment-length data, but this does not outperform the simpler power law model. For the earthquake frequency-moment data the gamma distribution provides the best overall fit to the data, implying a finite correlation length and a system near but below the critical point. In contrast, data from soft gamma ray repeaters show that the power law is the best fit, implying infinite correlation length and a system that is precisely critical. For the fault break data a significant break of slope is found instead at characteristic scale of 44 km, implying a typical seismogenic thickness of up to 22 km or so in west central Nevada. The exponent changes from 1.5 to -2.1, too large to be accounted for by changes in sampling for an ideal, isotropic fractal set. ¿ 2001 American Geophysical Union