1
$\begingroup$

I need to implement the following in python:

For a given discrete time series $Z_t$ ($t={0...T}$), find the smallest $t$ such that:

$$c\sum_{s=0}^t e^{[k(Z_t-Z_s)+m(t-s)]} \geq \frac{p^*}{1-p^*} $$

where $c,k,m$ are constants and

$$P^*=\int_0^{0.5} \frac{(1-2y)e^{-d/y}}{(1-y)^{1+d}y^{1-d}}dy=\int_{0.5}^{p^*} \frac{(2y-1)e^{-d/y}}{(1-y)^{1+d}y^{1-d}}dy $$

(d is another constant)

I need to implement this for each row of an $N\times w$ array $Z$.

I implemented the sum in equation 1 above utilizing the fact

$$ e^{[k(Zt-Zs)+m(t-s)]} = \frac{e^{kZt+mt}}{e^{kZs+ms}} $$

hence

$$\sum_{s=0}^t e^{[k(Z_t-Z_s)+m(t-s)]} =e^{kZ_t+m_t}*\sum_{s=0}^t\frac{1}{e^{kZ_s+m_s}} $$ which is reflected in use of cumsum in code below

def f(z):
    return ((1-2*z)*np.exp(-d/z))/(((1-z)**(1+d))*(z**(1-d)))
lhs=integrate.quad(f,0,0.5)

def rhs(p):
    return integrate.quad(-f,0.5, p)
p_star= fsolve(rhs-lhs,0.75) # will depend on time series only indrectly when d will be optimized

for i in np.arange(N):
    z=Z[i,:]  
    main = np.exp (kz+m*np.arange(w)) 
    cumsum_t=np.cumsum(1/main)
    final_sum= main*cumsum_t 
    t_solution= # index i where final_sum[i]> p_star/(1-p_star)) # not implemented yet

Is there any way I can vectorize this for $Z$? In this case $N$ is ~400,000 so vectorization would really help. f(z) will be fine with vector inputs but I don't think rhs(p) will be fine as it uses integration.

$\endgroup$
6
  • $\begingroup$ For the given sum $\sum_{s=0}^{t}$, what does it mean if $t \not \in \mathbb{Z}$ ? $\endgroup$
    – tqviet
    Oct 27, 2016 at 5:16
  • $\begingroup$ well, it is a discrete time series. $\endgroup$
    – dayum
    Oct 27, 2016 at 5:18
  • $\begingroup$ What are typical values of $p*$ and $d$? How accurate do you need $p*$ to be? If you want to vectorize the calculation of $p*$, I think the only way you can do it is by approximating the solution of the equation defining $p*$... $\endgroup$
    – GertVdE
    Oct 27, 2016 at 8:10
  • $\begingroup$ quad is incredibly slow even vectorizing this won't help much. Try a different integration method - I will dig up an old post I had that is way faster $\endgroup$
    – Matt
    Jun 17, 2017 at 16:12
  • $\begingroup$ stackoverflow.com/questions/37367688/… $\endgroup$
    – Matt
    Jun 17, 2017 at 16:15

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy