data21.rb

require 'statsample'



def f1(n)

	n = n * 3 + 1 if (n.odd?)
	return n

end


def f2(n)

	n = (n * 3 + 1) / 2 while (n.odd?)
	n /= 2 while (n.even?)
	return n

end

def stat(l)
	l = [0] if (l.empty?)
        t = t2 = 0
        l.each \
        {
		|x|
		t += x
		t2 += x ** 2
        }
        c = l.size
        a = t.to_f / c
        z = t2.to_f / c - a ** 2
        sd = Math.sqrt(z < 0 ? 0 : z)

	return a, sd, l.max.to_f, l.min.to_f
end

def stat2(l, t, n)
	return Hash[[["a#{n}", "sd#{n}", "mx#{n}"], stat(l)[0..2].map { |x| x / t }].transpose]
end

def d(s)
        c = s.split('').select { |x| x == '1' }.size
	d = c.to_f / s.length
        return d
        
end

def data(n)
	ns = n.to_s(2)
	nl = ns.length
	m = nl / 2
	
	nsh = ns[0..m]
	nsl = ns[m..-1]
	
	asdm1 = stat2(ns.split(/0+/).map { |x| x.length }, nl, 1)
	
	l1 = ns.split(/1+/)
	
	l1.shift
	asdm0 = stat2(l1.map { |x| x.length }, nl, 0)
	
	return {'d' => d(ns), 'dh' => d(nsh), 'dl' => d(nsl)}.merge(asdm1).merge(asdm0)
end

def dense(x, f)
	
	
	a = (0..99).to_a
	l = []
	x.times { l << a.delete_at(rand(a.size)) }
	s = '0' * 100
	l.each { |x| s[x, 1] = '1' }
	n = ('1' + s + '1').to_i(2)
	n2 = f.call(n)
	
	x = {}
	x['x1'] = n2.to_s(2).length.to_f / n.to_s(2).length
	x['x2'] = 1.0 / x['x1']
	
	x.merge!(data(n))
	x.merge!(Hash[data(n2).to_a.map { |k, v| ["#{k}_2", v] }])
	return x
end

def out(fn, a)
	return if (a.nil?)
	f = File.open(fn, 'a')
	f.puts(a.keys.join("\t")) if (f.size == 0)
	f.puts(a.values.join("\t"))
	f.close
	return a
end

def fit(l, y, lx)

	a = {}

	(lx + [y]).each { |x| a[x] = l.map { |b| b[x] }.to_vector() }

	ds = a.to_dataset()

	r = Statsample::Regression.multiple(ds, y)
#	$stderr.puts(r.summary)

	return r.coeffs.merge({'c' => r.constant})

end

def predict(l, y, z)
	l.each \
	{
		|x|
		t = z['c']
		
		(z.keys - ['c']).each { |k| t += z[k] * x[k] }
		x["#{y}_p"] = t
	}
end

def solve(l, y, x)
	z = fit(l, y, x)
#	$stderr.puts(z.inspect)
	predict(l, y, z)
end

def sum(l)

	t = 0.0
	l.each { |x| t += x }
	return t
end

def av(l)
	return nil if (l.empty?)
	return sum(l) / l.size
end

def corr(l, y1, yp)
	
	xav = av(l.map { |x| x[y1] })
	yav = av(l.map { |x| x[yp] })

	tx = ty = txy = e = 0.0
	l.each \
	{
		|z|
		
		x = z[y1]
		y = z[yp]
		
		txy += (x - xav) * (y - yav)
		tx += (x - xav) ** 2
		ty += (y - yav) ** 2
		e += (x - y) ** 2
	}
	r = txy / (Math.sqrt(tx) * Math.sqrt(ty))
	e /= l.size
	return r, e
	
end

def coef(l, y0, x0)
	r = nil
	begin
		solve(l, y0, x0)
		r, e = corr(l, y0, "#{y0}_p")
	rescue Statsample::Regression::LinearDependency
	end

	return {y0 => r}
end

def linear(fn)
	l = []
	100.times { |x| l << dense(x, method(fn)) }

	x0 = ['d', 'dh', 'dl', 'a0', 'sd0', 'mx0', 'a1', 'sd1', 'mx1']
	
	['x1', 'x2'].each \
	{
		|y|
		p(coef(l, y, x0))
		p(coef(l, y, x0.map { |x| "#{x}_2"}))
	}
	
	x0.each \
	{
		|y|
		
		p(coef(l, "#{y}_2", x0))
	}
	puts
end


linear(:f2)
linear(:f1)