You could use numpy.argmax:
ys = density(np.arange(9))
bb = np.argmax(ys)
aa = ys[bb]
This would compute the same values of aa
and bb
as the code you posted. However, this only finds the max among integer values of x
. If you look at Justin Peel's graph you'll see that the peak density can occur at some non-integer x-value
. So, to find a closer approximation to the peak density, use
xs = np.linspace(0,8,200)
ys = density(xs)
index = np.argmax(ys)
max_y = ys[index]
max_x = xs[index]