Comments (7)
from pacal.
@jszymon thanks for the quick response.
I used the version for python 3.
Please see my code below:
import time
from bisect import bisect
from scipy.optimize import minimize
grid = np.linspace(0, 7.61659217e-01,num=100)
def loglikelihood(param):
print("param:",param)
start_time = time.time()
#params = [sigma_bval, mu_tc, sigma_tc]
X0 = NormalDistr(0, param[0])
X1 = NormalDistr(param[1], param[2])
X2 = exp(X1)
X3 = X0 + X2
Y = [None]*len(grid)
for i in range(0, len(grid)):
if (i==0):
Y[i] = X3
else:
Y[i] = X3 + NormalDistr(0, grid[i])
def findYpdf(x, tsigt):
index = bisect(grid, tsigt)
pdf_l = Y[index-1].pdf(x)
pdf_u = Y[index].pdf(x)
#interpolate
pdf_final = pdf_l + (pdf_u - pdf_l)*(tsigt - grid[index-1])/(grid[index] - grid[index-1])
return pdf_final
LL = 0
for i in range(0, len(data)):
#start_time = time.time()
LL = LL - np.log(findYpdf(data.ix[i,'TC'], data.ix[i,'tsigt']))
#print("--- data",i,"is processed in %s seconds ---" % (time.time() - start_time))
print("--- processed in %s seconds ---" % (time.time() - start_time))
print(LL)
return LL
res = minimize(loglikelihood, [2.0,0,0.2], method='Nelder-Mead', tol=1e-6)
from pacal.
@jszymon I actually modified my code as follows, but the same error still raise.
import time
from bisect import bisect
from scipy.optimize import minimize
grid = np.linspace(0, 7.61659217e-01,num=100)
def loglikelihood(param):
print("param:",param)
start_time = time.time()
#params = [sigma_bval, mu_tc, sigma_tc]
X0 = NormalDistr(0, param[0])
X1 = NormalDistr(param[1], param[2])
X2 = exp(X1)
X3 = X0 + X2
Y = [None]*len(grid)
for i in range(0, len(grid)):
if (i==0):
Y[i] = X3
else:
Y[i] = X3 + NormalDistr(0, grid[i])
def findYpdf(x, tsigt):
try:
index = grid.tolist().index(tsigt)
pdf_final = Y[index].pdf(x)
except:
index = bisect(grid, tsigt)
pdf_l = Y[index-1].pdf(x)
pdf_u = Y[index].pdf(x)
#interpolate
pdf_final = pdf_l + (pdf_u - pdf_l)*(tsigt - grid[index-1])/(grid[index] - grid[index-1])
return pdf_final
LL = 0
for i in range(0, len(data)):
#start_time = time.time()
LL = LL - np.log(findYpdf(data.ix[i,'TC'], data.ix[i,'tsigt']))
#print("--- data",i,"is processed in %s seconds ---" % (time.time() - start_time))
print("--- processed in %s seconds ---" % (time.time() - start_time))
print(LL)
return LL
res = minimize(loglikelihood, [2.0,0,0.2], method='Nelder-Mead', tol=1e-6)
from pacal.
from pacal.
@jszymon
please see below:
import time
from bisect import bisect
params= [ 0.78440987, -0.01229896, -0.06878736]
grid = np.linspace(0, 7.61659217e-01,num=100)
tsigt = 0.10695393814651923
try:
index = grid.tolist().index(tsigt)
except:
index = bisect(grid, tsigt)
X0 = NormalDistr(0, params[0])
X1 = NormalDistr(params[1], params[2])
X2 = exp(X1)
X3 = X0 + X2
Y = [None]*len(grid)
for i in range(0, len(grid)):
if (i==0):
Y[i] = X3
else:
Y[i] = X3 + NormalDistr(0, grid[i])
pdf_l = Y[index-1].pdf(0.02049999999999841)
pdf_l
from pacal.
from pacal.
@jszymon good catch!
from pacal.
Related Issues (11)
- Assertion error when adding two r.vs (Normal) HOT 2
- Inverse CDF (PPF) fails closer to zero (outputs nan) HOT 2
- Messages about missing compiled routines should go to stderr HOT 2
- Requirements not specified correctly
- An import "numpy.fft.fftpack" is missing HOT 3
- On computing Pr( Y-const1 < X <= Y+const2) for random vars X and Y HOT 10
- Improving `Distr.hist()` HOT 3
- Log logistic function HOT 6
- Making a Fourier Transfers Program.
- DiracSegment & Segment with 0 probability ... can we skip them? HOT 1
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from pacal.