モジュール:Harmonic entropy

このモジュールについての説明文ページを モジュール:Harmonic entropy/doc に作成できます

local limits = require("Module:Limits")
local rat = require("Module:Rational")
local p = {}

-- Compute Harmonic Shannon entropy for an interval of `c` cents
-- `c`, `deviation`: in cents
-- `ratios`: an array of rational numbers
-- `norm`: a function of rational numbers
function p.harmonic_entropy(c, ratios, deviation, norm)
	norm = norm or function(ratio)
		return math.sqrt(rat.benedetti_height(ratio))
	end
	deviation = deviation or 1200 * math.log(1.01, 2)
	ratios = ratios
		or limits.integer_limit(200, function(ratio)
			if math.abs(rat.cents(ratio) - c) > 3 * deviation then
				return 1 / 0
			end
			return norm(ratio)
		end, 100)

	local function gaussian(x)
		return math.exp(-x * x / (2 * deviation * deviation)) / (deviation * math.sqrt(2 * math.pi))
	end

	local function weighted_gaussian(ratio)
		return gaussian(rat.cents(ratio) - c) / norm(ratio)
	end

	local q_norm = 0
	for _, ratio in pairs(ratios) do
		q_norm = q_norm + weighted_gaussian(ratio)
	end

	local function probability(ratio)
		return weighted_gaussian(ratio) / q_norm
	end

	local entropy = 0
	for _, ratio in pairs(ratios) do
		local p_i = probability(ratio)
		if p_i > 1e-5 then
			entropy = entropy - p_i * math.log(p_i)
		end
	end
	return entropy
end

return p