Este é meu código:
### BIBLIOTECAS
import scipy.special as sps
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import scipy as sp
from scipy.stats import norm
from scipy.stats import gamma
from math import exp
import operator
import csv
########################################
############### GERAR VALORES DA DISTRIBUIO GE ################
def rge(n, alpha, beta):
u = np.random.uniform(low=0, high=1, size=n)
x = -beta * np.log((-1 / alpha) * np.log(u))
return (x)
teste = rge(100, 1, 3)
plt.hist(teste)
plt.show()
temp = pd.DataFrame(teste)
########################################
### FUNCAO DE DISTRIBUICAO
# t[0] = alfa
# t[1] = beta
def fx(x, t):
prod = 1.0
for i in range(0, len(x)):
prod *= (t[0] / t[1]) * exp(- (x[i] / t[1])) * exp(-t[0] * exp(-(x[i] / t[1])))
return prod
##########################################
##### FUNCAO DE VEROSSIMILHANCA
def L(x, t):
return fx(x, t)
##########################################
### FUNCAO MCMC
def mcmc(N, k={"t0": 1, "t1": 1}, x=[]):
chute = {"t0": [1], "t1": [1]}
M = chute
hiper = {"t0": [0.1, 0.1], "t1": [0.1, 0.1]} # VALORES DOS HIPERPARAMETROS
contador = {"t0": [], "t1": []} # CONTADOR PARA TAXA DE ACEITACAO
thetas = M.keys()
for i in range(N - 1):
for j in thetas:
if j == "t0":
M[j].append(np.random.gamma(shape=M[j][-1], scale=k[j]))
lista = [[M[l][-1] for l in thetas], [M[l][-1] if l != j else M[l][-2] for l in thetas]]
t1 = gamma.pdf(M[j][-1], a=hiper[j][0], scale=hiper[j][1]) * L(x, lista[0]) * gamma.pdf(M[j][-2],
a=M[j][-1],
scale=k[j])
t2 = gamma.pdf(M[j][-2], a=hiper[j][0], scale=hiper[j][1]) * L(x, lista[1]) * gamma.pdf(M[j][-1],
loc=M[j][-2],
scale=k[j])
teste = (t1 / t2)
else:
M[j].append(np.random.gamma(shape=M[j][-1], scale=k[j]))
lista = [[M[l][-1] for l in thetas], [M[l][-1] if l != j else M[l][-2] for l in thetas]]
t1 = gamma.pdf(M[j][-1], a=hiper[j][0], scale=hiper[j][1]) * L(x, lista[0]) * gamma.pdf(M[j][-2],
a=M[j][-1],
scale=k[j])
t2 = gamma.pdf(M[j][-2], a=hiper[j][0], scale=hiper[j][1]) * L(x, lista[1]) * gamma.pdf(M[j][-1],
a=M[j][-2],
scale=k[j])
teste = (t1 / t2)
if (min(1, teste) < np.random.uniform(low=0, high=1)) or (np.isinf(teste)) or (np.isnan(teste)):
M[j][-1] = M[j][-2]
contador[j].append(0)
else:
contador[j].append(1)
print("Tamanho do theta 0 : " + str(len(M["t0"])))
print("\nTamanho do theta 1 : " + str(len(M["t1"])))
M = pd.DataFrame.from_dict(M)
contador = pd.DataFrame.from_dict(contador)
cont = contador.apply(sum)
print(cont)
return (M)
N = int(input("Entre com o N: "))
MP = mcmc(N=N, x=temp)
print(MP)
E está gerando o seguinte erro:
Traceback (most recent call last):
File "/home/karlla/Documentos/Documentos/Mestrado/2° semestre/Estatistica computacional /Trabalho-artigo2/gerar.py", line 116, in <module>
MP = mcmc(N=N, x=temp)
File "/home/karlla/Documentos/Documentos/Mestrado/2° semestre/Estatistica computacional /Trabalho-artigo2/gerar.py", line 74, in mcmc
t1 = gamma.pdf(M[j][-1], a=hiper[j][0], scale=hiper[j][1]) * L(x, lista[0]) * gamma.pdf(M[j][-2],
File "/home/karlla/Documentos/Documentos/Mestrado/2° semestre/Estatistica computacional /Trabalho-artigo2/gerar.py", line 53, in L
return fx(x, t)
File "/home/karlla/Documentos/Documentos/Mestrado/2° semestre/Estatistica computacional /Trabalho-artigo2/gerar.py", line 44, in fx
prod *= (t[0] / t[1]) * exp(- (x[i] / t[1])) * exp(-t[0] * exp(-(x[i] / t[1])))
File "/home/karlla/anaconda2/lib/python2.7/site-packages/pandas/core/series.py", line 78, in wrapper
"{0}".format(str(converter)))
TypeError: cannot convert the series to <type 'float'>
Não estou conseguindo arrumar isso.