All pastes #433347 Raw Edit

aegray

public text v1 · immutable
·published 2007-04-10 19:03 UTC
rendered paste body




###Generate list of x

x = NULL
n = 50
pixprime = NULL
pix = NULL
pixprimeratio = NULL
pixratio = NULL
accepted = 0
acceptanceratio = 0
mag = NULL


for(i in 1:n)
{

spinran = runif(1,0,1)


if(spinran > .5)
{
spinstate = 1
}

if(spinran <= .5)
{
spinstate = -1
}

x<-c(x, spinstate)

}



###Burn for 10,000 times and run for 100,000

T = 100000

for(i in 1:T)
{

	###Change Random Spin




	index = ceiling(runif(1,0,50))

	xprime = x

	xprime[index] = xprime[index]*-1


	###Computing the target distribution and metropolis ratio
	pix <- NULL
	pixprime <- NULL
	for(i in 1:n-1)
	{
		pix = c(pix, 2*x[i]*x[i+1])
		pixprime = c(pixprime, 2*xprime[i]*xprime[i+1])

	}

	pixratio = exp(sum(pix))
	pixprimeratio = exp(sum(pixprime))



	MetRatio <- pixprimeratio/pixratio


	###test answer

	U <- runif(1,0,1)
	if (U <= MetRatio)
	{
		x = xprime
		if(T>10000)
		{
			accepted = accepted + 1
		}
	}

	

	if(T>10000)
	{
	mag = c(mag, sum(x))
	}




}
acceptanceratio = accepted/90000
acceptanceratio

acf(mag)
hist(mag)