From: Chris Pine <nemo@...>
Date: 2003-04-19T05:33:39+09:00
Subject: Re: Biased weighted random?

All right...

It's not pretty, but it works and you only have to do it once, right at the
beginning of your endless random shuffle thing.

The idea is to use the "bonehead" algorithm, but not with your target
probabilities as the weights, since that won't work.

Instead, use this program.  Enter in your weights into the `weights' array
(inside the call to `normalize' and run it.

Take the new weights to use (the first bunch of numbers) and put them into
your "bonehead" algorithm as the weights (assuming the probabilities shown
(second bunch of numbers) is close enough for you... if not, change the `50'
in `50.times' to something bigger).

Should work much of the time...

Chris


-----------------------------


require 'matrix'


def normalize vector
  d = vector.inject(0.0){|sum,x| sum+x}  #  common denominator
  vector.map{|x| x/d}
end

def realProbVector vector
  n = vector.length

  #  Matrix is an array of cols.
  cols = Array.new(n){Array.new(n)}
  for x in 0...n
    for y in 0...n
      cols[x][y] = 0
      next if x == y
      cols[x][y] = vector[x]/(1-vector[y])
    end
  end

  transitionMatrix = Matrix.columns(cols)

  newMatrix = transitionMatrix

  5.times {newMatrix = newMatrix * newMatrix}

  newMatrix.row(0)
end


targets = normalize([2,1,1])
weights = targets.dup
probs   = realProbVector(weights)

50.times do
  weights.each_with_index do |x,idx|
    weights[idx] += targets[idx] - probs[idx]
  end
  weights = normalize(weights)
  probs = realProbVector(weights)
end


weights.to_a.each do |prob|
  puts prob.to_s
end
puts
probs.to_a.each do |prob|
  puts prob.to_s
end