Oh you can absolutely average rotation quaternions, even weighted. You just use regular online computation of mean, express it in terms of linear mix operations (x k + (1-k) y) and then use quaternion slerp for the linear mix:

  mean = zero_rotation_quaternion
  wsum = 0
  for w,qt in zip(weights, quaternions):
    wsum += w
    # regular online mean computation:
    # delta = qt - mean
    # mean := mean + delta*w/wsum
    # express as linear mix:
    # mean := mean + qt*w/wsum - mean*w/wsum
    # mean := mean * (1-w/wsum) + qt*w/wsum
    # i.e. mean is linear mix between mean and qt, with k=w/sum
    mean = slerp(mean, qt, w/wsum)
  return mean