Rich: one approach is to fit the data to a gamma distribution, and take the mean and variance from that. I've found Bolker's book, Ecological Models & Data in R, to be helpful in understanding the gamma distribution. If you're not familiar with Newman et al. Estimating mean and variance for environmental studies with below detection limit observations Water Resources Bulletin 25, pp. 905-916, it is mandatory reading. Below is some R code provided by a statistician, Jim Baldwin, for estimating [NO3] mean and variance in an experiment involving three tree species and two levels of plant diversity. NO3 concentrations have many of the properties you describe for metal concentrations. The code has provisions for dealing with left-censored data, but the exact details escape me at the moment. --Seth

# Estimating gamma distribution parameters given the data is
# left-censored at 0.04.

# Rather than estimate the shape and scale parameters and then
# determining the mean and variance, the mean and variance are
# found directly by re-parametrizing the gamma density function.

# Get original data

  setwd("c:\\documents and settings\\jbaldwin\\desktop\\seth")

  g = read.csv("Task1404_Treats_3-4 NO3 min set.csv")

# Define function returning the log of the likelihood

logL = function(z,x) {

  shape = z[1]^2/z[2] # Shape parameter from mean (z[1]) and variance (z[2])

  scale = z[2]/z[1]   # Scale parameter from mean (z[1]) and variance (z[2])

  # Non-censored data values

  xx = x[(!is.na(x)) & (x > 0.04)]

  # Number of censored data values

  n = length(x[(!is.na(x)) & (x <= 0.04)])

  # Return log of the likelihood

  n*pgamma(0.04,shape=shape,scale=scale,log.p=TRUE) +

    sum(dgamma(xx,shape=shape,scale=scale,log=TRUE))
}

# MLE function for mean and variance
mle.censored = function(species,treatment) {

  # Extract data for specified species and treatment
  x = g$NO3[g$KEYSPECI==species & g$TREATMNT==treatment]

  # Get maximum likelihood estimates
  mle.g = optim(c(1,1),logL,method = "L-BFGS-B",lower=0.0001,upper=Inf,
    hessian=TRUE,x=x,control=list(fnscale=-1))

  # Calculate estimate of covariance matrix
  mle.g$VarCovMatrix = -solve(mle.g$hessian)

  # Save a copy of the data
  mle.g$x = x

  # Return the data frame with all of the results
  mle.g

}

# Get estimates for each combination of species and treatment
  mle13 = mle.censored(1,3)
  mle14 = mle.censored(1,4)
  mle23 = mle.censored(2,3)
  mle24 = mle.censored(2,4)
  mle33 = mle.censored(3,3)
  mle34 = mle.censored(3,4)

# Collect all terms
  # Maximum likelihood estimates of population means
  mu.mle = c(mle13$par[1],mle14$par[1],mle23$par[1],mle24$par[1],
             mle33$par[1],mle34$par[1])
  # Maximum likelihood estimates of population variances
  var.mle = c(mle13$par[2],mle14$par[2],mle23$par[2],mle24$par[2],
              mle33$par[2],mle34$par[2])
  # Estimate of standard error of estimates of population means
  mu.se = c(mle13$VarCov[1,1],mle14$VarCov[1,1],mle23$VarCov[1,1],
            mle24$VarCov[1,1],mle33$VarCov[1,1],mle34$VarCov[1,1])^0.5
  # Estimate of standard error of estimates of population variances
  var.se = c(mle13$VarCov[2,2],mle14$VarCov[2,2],mle23$VarCov[2,2],
            mle24$VarCov[2,2],mle33$VarCov[2,2],mle34$VarCov[2,2])^0.5
  # Species and treatment id
  id = c("s1 t3","s1 t4","s2 t3","s2 t4","s3 t3","s3 t4")

# Plot means and approximate 95% confidence intervals
  plot(c(min(mu-1.96*mu.se),max(mu+1.96*mu.se)),c(0.5,6.5),type="n",
    xlab="Mean NO3 concentration",ylab="",axes=FALSE)
  axis(1)
  axis(2,c(1:6),id,las=1)
  box()
  for (i in c(1:6)) {
      points(mu[i],i,pch=16)
      lines(c(mu[i]-1.96*mu.se[i],mu[i]+1.96*mu.se[i]),c(i,i))
  }




This electronic message contains information generated by the USDA solely for the intended recipients. Any unauthorized interception of this message or the use or disclosure of the information it contains may violate the law and subject the violator to civil or criminal penalties. If you believe you have received this message in error, please notify the sender and delete the email immediately.

	[[alternative HTML version deleted]]

