Улучшить скорость работы петли
Я пытаюсь создать пример для кода Монте-Карло моей цепочки Маркова, используя pyMC. Таким образом, с выбранными параметрами модели, каждый раз, когда вывод создается путем вызова getLensing Экземпляр из класса NFW и по сравнению с наблюдаемыми данными. Моя проблема в том, что мой код очень медленный, когда он вычисляет параметр модели. У меня есть, например, 24000 точек данных, а затем для каждого из них у меня есть распределение вероятностей, например obj_pdf- который я маргинализировать над ним (интегрировать) во внутреннем цикле. поэтому каждый раз, по крайней мере, требуется час, чтобы вычислить все выходные данные модели.
import numpy as np
z=np.arange(0,1.5,0.001)
z_h=0.15
for j in range(pos.shape[0]):
value1=0;value2=0
pdf=obj_pdf[j,:]/sum(obj_pdf[j,:])
for i in range(len(z)):
if (z[i]>z_h) :
g1,g2=nfw.getLensing( pos[j,:], z[i])
value1+=g1*pdf[i]
value2+=g2*pdf[i]
if (j<1):
value=np.array([value1,value2])
else:
value=np.vstack((value, np.array([value1,value2])))
поэтому, если я захочу повторно сэмплировать входные параметры, например, 100000 раз, у меня будут месяцы для расчета MCMC. Есть ли какой-нибудь умный способ ускорить мой код и циклы? Нужно ли использовать что-то вроде numpy.vectorize или все же это не улучшит скорость моего кода? Как насчет cythonэто увеличит производительность кода? Если это помогает, как это работает?
Я побежал python -m cProfile mycode.py чтобы увидеть, что было вызвано тем, что мой код работает медленно, и это были результаты:
12071 0.004 0.000 0.004 0.000 {min}
2 0.000 0.000 0.000 0.000 {next}
1 0.000 0.000 0.000 0.000 {numexpr.interpreter._set_num_threads}
8 0.002 0.000 0.002 0.000 {numpy.core.multiarray.arange}
132424695 312.210 0.000 312.210 0.000 {numpy.core.multiarray.array}
73498 3.933 0.000 3.933 0.000 {numpy.core.multiarray.concatenate}
99151506 201.497 0.000 201.497 0.000 {numpy.core.multiarray.copyto}
99151500 164.303 0.000 164.303 0.000 {numpy.core.multiarray.empty_like}
28 0.000 0.000 0.000 0.000 {numpy.core.multiarray.empty}
2 0.000 0.000 0.000 0.000 {numpy.core.multiarray.set_string_function}
1 0.000 0.000 0.000 0.000 {numpy.core.multiarray.set_typeDict}
1 0.000 0.000 0.000 0.000 {numpy.core.multiarray.where}
14 0.000 0.000 0.000 0.000 {numpy.core.multiarray.zeros}
14 0.000 0.000 0.000 0.000 {numpy.core.umath.geterrobj}
7 0.000 0.000 0.000 0.000 {numpy.core.umath.seterrobj}
270 0.000 0.000 0.000 0.000 {numpy.lib._compiled_base.add_docstring}
6 0.011 0.002 0.011 0.002 {open}
1 0.000 0.000 0.000 0.000 {operator.div}
2 0.000 0.000 0.000 0.000 {operator.mul}
1918 0.000 0.000 0.000 0.000 {ord}
2 0.000 0.000 0.000 0.000 {posix.WEXITSTATUS}
2 0.000 0.000 0.000 0.000 {posix.WIFEXITED}
1 0.000 0.000 0.000 0.000 {posix.WIFSIGNALED}
9 0.002 0.000 0.002 0.000 {posix.access}
3 0.000 0.000 0.000 0.000 {posix.close}
5 0.002 0.000 0.002 0.000 {posix.fdopen}
1 0.002 0.002 0.002 0.002 {posix.fork}
4 0.000 0.000 0.000 0.000 {posix.getcwd}
6 0.000 0.000 0.000 0.000 {posix.getpid}
1 0.000 0.000 0.000 0.000 {posix.getuid}
1 0.000 0.000 0.000 0.000 {posix.listdir}
6 0.000 0.000 0.000 0.000 {posix.lstat}
4 0.043 0.011 0.043 0.011 {posix.open}
2 0.000 0.000 0.000 0.000 {posix.pipe}
2 0.004 0.002 0.004 0.002 {posix.popen}
1 0.007 0.007 0.007 0.007 {posix.read}
205 0.059 0.000 0.059 0.000 {posix.stat}
3 0.000 0.000 0.000 0.000 {posix.sysconf}
2 0.000 0.000 0.000 0.000 {posix.uname}
4 0.004 0.001 0.004 0.001 {posix.unlink}
3 0.000 0.000 0.000 0.000 {posix.urandom}
1 0.000 0.000 0.000 0.000 {posix.waitpid}
1 0.000 0.000 0.000 0.000 {pow}
1522 0.004 0.000 0.004 0.000 {range}
73 0.000 0.000 0.000 0.000 {repr}
99151501 2102.879 0.000 6380.906 0.000 {scipy.integrate._quadpack._qagse}
1776 0.002 0.000 0.002 0.000 {setattr}
32 0.000 0.000 0.000 0.000 {sorted}
24500 18.861 0.001 18.861 0.001 {sum}
184 0.000 0.000 0.000 0.000 {sys._getframe}
1 0.000 0.000 0.000 0.000 {sys.getfilesystemencoding}
2 0.000 0.000 0.000 0.000 {sys.settrace}
1 0.000 0.000 0.000 0.000 {tables.utilsextension._broken_hdf5_long_double}
1 0.000 0.000 0.000 0.000 {tables.utilsextension.blosc_compressor_list}
2 0.000 0.000 0.000 0.000 {tables.utilsextension.get_hdf5_version}
1 0.000 0.000 0.000 0.000 {tables.utilsextension.get_pytables_version}
2 0.000 0.000 0.000 0.000 {tables.utilsextension.which_lib_version}
27 0.000 0.000 0.000 0.000 {thread.allocate_lock}
6 0.000 0.000 0.000 0.000 {thread.get_ident}
4 0.000 0.000 0.000 0.000 {thread.start_new_thread}
1 0.000 0.000 0.000 0.000 {time.localtime}
2 0.000 0.000 0.000 0.000 {time.time}
105 0.000 0.000 0.000 0.000 {unichr}
229 0.000 0.000 0.000 0.000 {vars}
49300 2.127 0.000 2.127 0.000 {zip}
1 ответ
Вот код Я был бы удивлен, если бы время шло от 60 до 59 минут.
import numpy as np
z_h=0.15
z=np.arange(z_h, 1.5,0.001) #start the range from what you need (not exactly
z=z[1:] # needed because you said if (z[i]>z_h), range gives (z[i]>=z_h)
value=np.array([])
for j in range(pos.shape[0]):
value1=0;value2=0
pdf=obj_pdf[j,:]/sum(obj_pdf[j,:])
posj=pos[j,:] #precalculate
for i,zi in enumerate(z): #use enumerate if you need value and index
g1,g2=nfw.getLensing( posj, zi)
value1+=g1*pdf[i]
value2+=g2*pdf[i]
value=np.append(value, np.array([value1,value2])) # use a proper append function
Как и другие, я предполагаю, что getLensing пожирает ваши циклы процессора.
Согласно первому ответу, np.vectorize не будет ускорять вашу функцию.