How to compute wavelet coherence of two given timeseries in python?

Question

I want to compute wavelet coherence of two given timeseries in python. I tried working with pycwt but when calculating the cross-wavelet transform I find trouble due to shortness of data or too much of a trending (both things have already been taken care of so there is not much I can do). I compute the cross-wavelet manually and proceed to do the same with wavelet coherence, but I am in need of a smoothing function and it all gets messy and complicated. I don't think there is a simple way to do this, but maybe somebody have a clearer idea. I am following the paper "Partial wavelet coherence as a robust method for assesment of neurovascular coupling in neonates with hypoxic ischemic encephalopathy", by T. Hermans et al. I am also aware of a very useful R implementation of the wavelet coherence that computes it very easily, but I am asked to do the Python implementation exactly. Thank you everybody!

Tried pycwt package and more things on the way... I have been working on this for a while now. Code attached.

1st try, where TS is my timeseries:

import numpy as np
import pycwt as wv

#Parameters:
N = len(TS['signal1']) # Length of time series.
mother = Morlet(6)
dt = 1
s0 = float(2 * dt)
dj = 1/12
J = round(np.log2((N*dt)/s0) / dj)
#J = 7 / dj
scales = s0 * 2**(dj * np.arange(J))

cwt_y1 = wv.cwt(TS['signal1'], dt, dj, s0, J, wavelet='morlet', freqs=None)
cwt_y2 = wv.cwt(TS['signal2'], dt, dj, s0, J, wavelet='morlet', freqs=None)

# Compute cross-wavelet transform
cross_wvt, _, _ = wv.wct(TS['signal1'], TS['signal2'], dt, dj, s0, J, wavelet='morlet', freqs=None)

# Compute wavelet coherence
coherence, phases, _, _, _ = wv.wct(TS['signal1'], TS['signal2'], dt, dj, s0, J, wavelet='morlet', freqs=None, coherence=True)

# Normalize coherence matrix
coherence /= np.max(coherence)

Error: cross_wvt, _, _ = wv.wct(TS['ave_TOI1'], TS['TOTPOW08_BIS'], dt, dj, s0, J, wavelet='morlet', freqs=None) Warning: Cannot place an upperbound on the unbiased AR(1). Series is too short or trend is to large.

I was expecting a computation of wavelet coherence's matrix, as in R quicker implementation: wtc(signal1, signal2, dj=dj, J1=J, s0=s0, mother=mother, nrands=iterations)