вычисление функции взаимной корреляции в xarray

У меня есть набор данных res_1 с

Dimensions:    (space: 726, time: 579)
Coordinates:
  * space      (space) MultiIndex
  - latitude   (space) float64 -90.0 -82.5 -82.5 -82.5 -82.5 -82.5 -82.5 ...
  - longitude  (space) float64 0.0 0.0 60.0 120.0 180.0 240.0 300.0 0.0 30.0 ...
  * time       (time) datetime64[ns] 1980-06-01 1980-06-02 1980-06-03 ...
  Data variables:
       mx2t       (time, space) float64 -1.768 -0.6035 -1.286 -1.291 1.144 ...
       dayofyear  (time) int64 153 154 155 156 157 158 159 160 161 162 163 164 ...

пространственная переменная содержит пары широты и долготы. I для расчета функции взаимной корреляции

cij = (avg(mx2t(t-tau , i) * mx2t(t , j)) - avg(mx2t(t-tau , i))*avg(mx2t(t , j)))/(std(mx2t(t-tau , i))*std(mx2t , j) )

где avg - математическое ожидание, а std - стандартное отклонение, i и j проходят по всем элементам в пространственной координате, а tau изменяется от 0 до 200. для этого я определил функцию

def c_out(i) :
    c1=[]
    c = np.empty(726)
    c.fill(-2.0)
    c[i]=0.0
    for j in list(range(726)):
        if i != j :
            rdi = res_1.sel(space = coord[i]).to_dataframe()
            rdj = res_1.sel(space = coord[j]).to_dataframe()
            rdi['tj'] = rdj['t']
            for tau in list(range(200)):
                rdi['mx2t_stau'] = rdi['t'].shift(tau)
                rdf = rdi.dropna()
                rdf1 = rdf.loc[pd.date_range('1982-01-01' , '1982-12-31')]
                ctemp = ((rdf1['tj']*rdf1['mx2t_stau']).mean() - rdf1['tj'].mean() * rdf1['mx2t_stau'].mean()/(rdf1['tj'].std()*rdf1['mx2t_stau'].std())
                if ctemp > c[j] :
                     c[j] = ctemp
    return c    

и я использую joblib, чтобы вычислить его параллельно, используя

cij = Parallel(n_jobs=28 )(delayed(c_out)(i)for i in list(range(726))

Я хотел бы знать, есть ли простой или (/ и) более эффективный способ выполнения тех же вычислений в xarray?


person aniruddha    schedule 08.08.2017    source источник
comment
Что значит лучший способ сделать это? С точки зрения вычислительной производительности или алгоритмической реализации? Попробуйте пересмотреть свой вопрос, чтобы конкретизировать, в чем вам нужна помощь.   -  person jhamman    schedule 11.08.2017