Преобразование Фурье в Pandas

#python #pandas #dataframe #scipy #фурье

Есть DataFrame примерно такого содержания:

Time     Data
95114    1012.0
96078    1108.0
97033    1132.0
98015    1076.0
...


Состоит из 100 строк, где:


Time - время измерения в микросекундах
Data - измеряемое значение. 


Мне необходимо найти ДПФ этой зависимости. Все бы ничего, только частота измерения
непостоянна. Она немного плывет со временем случайным образом. Т.е. мои потуги в духе:

xf = fftfreq(len(df['Time']),1000) # 1000 взята как примерный период измерений
i = xf > 0
yf = fft(df['Data'])


Дают немного неожиданный результат. Я так понимаю, что дело в частоте. Как правильно
ее искать? Интерполировать свой дискретный ряд и потом уже находить его ДПФ? 

В итоге мне нужно будет найти площадь между определенными частотами спектра. 
    

Ответы

Ответ 1

Можно воспользоваться методом DataFrame.resample(), чтобы преобразовать DataFrame
во временной ряд с единым интервалом (в качестве единого интервала возьмем средний
интервал между замерами в микросекундах) и после этого интерполировать недостающие
значения.

Для этого надо конвертировать df.Time в datetime dtype и установить этот столбец
в качестве индекса:

In [62]: freq = int(df['Time'].diff().mean()) # средняя частота измерения в микросекундах

In [63]: freq
Out[63]: 967

In [64]: df = (df.set_index(pd.to_datetime(df.pop('Time'), unit='us'))
                 .resample('{}U'.format(freq))
                 .mean()
                 .interpolate())

In [65]: df
Out[65]:
                              Data
Time
1970-01-01 00:00:00.094766  1012.0
1970-01-01 00:00:00.095733  1108.0
1970-01-01 00:00:00.096700  1132.0
1970-01-01 00:00:00.097667  1076.0


После этого можно попробовать воспользоваться преобразованиями Фурье.

Что делает данная операция `img[52:308,52:308]`

#python #изображения #numpy #scipy #pil

Изучаю питон, и никак не пойму, что делают строки:

img=img[52:308,52:308] 
print img.shape


В коде:  

#Fit the data into a 256*256 format and save it as a .tif file
#@Sukriti Paul

from PIL import Image
from skimage import data, io, filters
import numpy as np
from scipy import ndimage
from matplotlib import pyplot as plt
from sklearn import manifold, datasets

n=10 #number of images
for i in range (1,n+1):
    path=str(i)+'.png'
    img= io.imread('/Users/sukritipaul/Desktop/D004/test_annotate/'+path)
      #Subsection of the image
    img=img[52:308,52:308]
    print img.shape
    path=str(i)+'.tif'
    io.imsave('/Users/sukritipaul/Desk

    

Ответы

Ответ 1

Это обрезка полей изображения (cropping) - все пиксели по X и Y обрезаются таким
образом - выбрать все пиксели от 52-го до 308-го по соответствующей оси и пропустить
все остальные.

Пример в Numpy:

In [267]: import numpy as np

In [268]: a = np.arange(36).reshape(6,6)

In [269]: a
Out[269]:
array([[ 0,  1,  2,  3,  4,  5],
       [ 6,  7,  8,  9, 10, 11],
       [12, 13, 14, 15, 16, 17],
       [18, 19, 20, 21, 22, 23],
       [24, 25, 26, 27, 28, 29],
       [30, 31, 32, 33, 34, 35]])


обрежем по одному пикселю:

In [270]: a = a[1:5, 1:5]

In [271]: a
Out[271]:
array([[ 7,  8,  9, 10],
       [13, 14, 15, 16],
       [19, 20, 21, 22],
       [25, 26, 27, 28]])


Фактически это обычная операция "среза" (slicing), которая в Numpy массивах может
быть одновременна выполнена по нескольким осям (измерениям).

Для N-мерной матрицы срезы можно делать сразу по всем измерениям:

arr[, , ..., ]




Атрибут numpy.ndarray.shape возвращает размерность массива в виде кортежа:

In [295]: a.shape
Out[295]: (4, 4)

In [296]: np.arange(48).reshape(8,6).shape
Out[296]: (8, 6)

In [297]: np.arange(2*3*4).reshape(2,3,4).shape
Out[297]: (2, 3, 4)


Ответ 2

Вырезают фрагмент кадра размером 256x256 из, видимо, бОльшего формата, начиная с
точки (52, 52)

Python. Как быстро посчитать расстояние между всеми точками в массивах?

#python #массивы #scipy #spatial #haversine

Есть 733  географических точек F (массив 2*733)
Есть 22 000 географических точек S (массив 2*22 000)

Как быстро посчитать расстояние от каждой точки из F до каждой точки S.

Этот код работает очень долго, около семи минут. Почему?

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda u, v: haversine(u, v))


А если вычислять евклидову метрику, то быстро.

import scipy.spatial
D2 = scipy.spatial.distance.cdist(F,S,'euclidean')


Вычислять расстояние функцией  тоже долго:

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda u,v:geopy.distance.vincenty(u,v).km) 


Но дело не в сложности функции вычисления расстояния, т.к даже если я буду считать
сумму разности координат вычисления займут много времени:

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda  u,v:(u-v).sum())

    

Ответы

Ответ 1

Вот рабочий пример.

Данные для примера: http://samplecsvs.s3.amazonaws.com/SacramentocrimeJanuary2006.csv

import numpy as np
import pandas as pd
from sklearn.neighbors import DistanceMetric

# http://scikit-learn.org/stable/modules/generated/sklearn.neighbors.DistanceMetric.html
dist = DistanceMetric.get_metric('haversine')

url = 'http://samplecsvs.s3.amazonaws.com/SacramentocrimeJanuary2006.csv'
df = pd.read_csv(url, usecols=['latitude','longitude'])


первые 5 записей:

In [213]: df.head()
Out[213]:
    latitude   longitude
0  38.550420 -121.391416
1  38.473501 -121.490186
2  38.657846 -121.462101
3  38.506774 -121.426951
4  38.637448 -121.384613


увеличим DataFrame в три раза, чтобы получить примерно такой же объем данных как
у вас в вопросе:

df = pd.concat([df] * 3, ignore_index=True).sample(frac=1)
# split our data set into two
F, S = np.split(df, [int(len(df)*0.033)])


получилось два DataFrame'а:

print('shape of [F]: {}'.format(F.shape))
print('shape of [S]: {}'.format(S.shape))


выдает:

shape of [F]: (750, 2)
shape of [S]: (22002, 2)


считаем расстояние:

earth_radius = 6371

D = dist.pairwise(np.radians(F), np.radians(S)) * earth_radius


PS haversine считается по формуле:

2 * arcsin(sqrt(sin^2(0.5*dx) * cos(x1) * cos(x2) * sin^2(0.5*dy)))


и возвращает результат в радианах, который надо умножить на радиус Земли (в км.),
чтобы получить расстояние в километрах

Результат:

In [208]: D
Out[208]:
array([[  3.86258192,  12.14112631,  16.1702106 , ...,  11.2399247 ,  13.71076021,
  2.59575523],
       [  3.19765341,  15.72580688,  20.35963524, ...,  15.38247654,  17.85975202,
  4.60552864],
       [  6.15251463,  19.94251384,  19.77488568, ...,  17.0637875 ,  19.15736473,
  6.55354701],
       ...,
       [  7.85715424,   9.52674314,  11.13189495, ...,   6.81544901,   9.16104377,
  5.58452969],
       [ 14.02227805,   5.08370143,   6.47846773, ...,   0.57188025,   3.04569744,
 11.81577747],
       [ 17.84804533,   3.93791815,   6.91865107, ...,   3.8743883 ,   2.94376778,
 15.72668196]])

In [209]: D.shape
Out[209]: (750, 22002)


замер скорости - около 6.5 секунд (для 16.5 миллионов пар точек) на моем стареньком
ноуте:

In [210]: %timeit dist.pairwise(np.radians(F), np.radians(S)) * earth_radius
6.45 s ± 156 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)


Проверка формулы:

In [204]: lyon = np.array([[45.7597, 4.8422]])

In [205]: paris = np.array([[48.8567, 2.3508]])

In [206]: dist.pairwise(np.radians(lyon), np.radians(paris))  * earth_radius
Out[206]: array([[ 392.21671781]])

Преобразование Фурье в Pandas

Есть DataFrame примерно такого содержания:
Time Data 95114 1012.0 96078 1108.0 97033 1132.0 98015 1076.0 ...
Состоит из 100 строк, где:
Time - время измерения в микросекундах Data - измеряемое значение.
Мне необходимо найти ДПФ этой зависимости. Все бы ничего, только частота измерения непостоянна. Она немного плывет со временем случайным образом. Т.е. мои потуги в духе:
xf = fftfreq(len(df['Time']),1000) # 1000 взята как примерный период измерений i = xf > 0 yf = fft(df['Data'])
Дают немного неожиданный результат. Я так понимаю, что дело в частоте. Как правильно ее искать? Интерполировать свой дискретный ряд и потом уже находить его ДПФ?
В итоге мне нужно будет найти площадь между определенными частотами спектра.

Ответ

Можно воспользоваться методом DataFrame.resample(), чтобы преобразовать DataFrame во временной ряд с единым интервалом (в качестве единого интервала возьмем средний интервал между замерами в микросекундах) и после этого интерполировать недостающие значения.
Для этого надо конвертировать df.Time в datetime dtype и установить этот столбец в качестве индекса:
In [62]: freq = int(df['Time'].diff().mean()) # средняя частота измерения в микросекундах
In [63]: freq Out[63]: 967
In [64]: df = (df.set_index(pd.to_datetime(df.pop('Time'), unit='us')) .resample('{}U'.format(freq)) .mean() .interpolate())
In [65]: df Out[65]: Data Time 1970-01-01 00:00:00.094766 1012.0 1970-01-01 00:00:00.095733 1108.0 1970-01-01 00:00:00.096700 1132.0 1970-01-01 00:00:00.097667 1076.0
После этого можно попробовать воспользоваться преобразованиями Фурье.