FilesExpand file tree

 master

/

algebra.py

Copy path

Blame

More file actions

Blame

More file actions

Latest commit

 

History

History

History

61 lines (45 loc) · 1.44 KB

 master

/

algebra.py

Copy path

Top

File metadata and controls

• Code
• Blame

61 lines (45 loc) · 1.44 KB

Raw

Copy raw file

Download raw file

Open symbols panel

Edit and raw actions

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

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

59

60

61

from matplotlib.finance import quotes_historical_yahoo

from datetime import date

import numpy as np

import matplotlib.pyplot as plt

from scipy import fftpack

from scipy import signal

from matplotlib.dates import DateFormatter, DayLocator, MonthLocator

from scipy import optimize

today = date.today()

start = (today.year - 1, today.month, today.day)

quotes = quotes_historical_yahoo("QQQ", start, today)

quotes = np.array(quotes)

dates = quotes.T[0]

qqq = quotes.T[4]

y = signal.detrend(qqq)

alldays = DayLocator()

months = MonthLocator()

month_formatter = DateFormatter("%b %Y")

fig = plt.figure()

ax = fig.add_subplot(211)

ax.xaxis.set_minor_locator(alldays)

ax.xaxis.set_major_locator(months)

ax.xaxis.set_major_formatter(month_formatter)

amps = np.abs(fftpack.fftshift(fftpack.rfft(y)))

amps[amps < amps.max()] = 0

def residuals(p, y, x):

A,k,theta = p

err = y-A * np.sin(2* np.pi* k * x + theta)

return err

filtered = -fftpack.irfft(fftpack.ifftshift(amps))

N = len(qqq)

f = np.linspace(-N/2, N/2, N)

p0 = [filtered.max(), f[amps.argmax()]/N, np.pi/3]

plsq = optimize.leastsq(residuals, p0, args=(filtered, dates))

p = plsq[0]

print p

plt.plot(dates, y, 'o', label="detrended")

plt.plot(dates, filtered, label="filtered")

plt.plot(dates, p[0] * np.sin(2 * np.pi * dates * p[1] + p[2]), '^', label="fit")

fig.autofmt_xdate()

plt.legend()

ax2 = fig.add_subplot(212)

plt.plot(f, amps, label="transformed")

plt.legend()

plt.show()