Skip to content

Commit 6a01691

Browse files
author
gigi
committed
update
1 parent e65c50b commit 6a01691

149 files changed

Lines changed: 3424 additions & 0 deletions

File tree

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

NumPy/Code/ch10code/README

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
NumPy Beginner's Guide
2+
3+
Chapter 10 When Numpy is not enough: Scipy and beyond
4+
5+
The example code requires SciPy and Matplotlib.
6+
7+
Run the code examples with
8+
9+
python <name of script>
10+
11+
frequencies.py Demonstrates FFT transforms.
12+
gaussquad.py Demonstrates numerical integration.
13+
images.py Demonstrates image manipulation.
14+
optfit.py Demonstrates optimization.
15+
pair.py Demonstrates two samples statistical tests.
16+
scipyio.py Demonstrates Matlab/Octave interoperability.
17+
sincinterp.py Demonstrates interpolation.
18+
statistics.py Demonstrates statistical functionality.
19+
trend.py Demonstrates the detrend function.

NumPy/Code/ch10code/a.mat

248 Bytes
Binary file not shown.

NumPy/Code/ch10code/algebra.py

Lines changed: 61 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,61 @@
1+
from matplotlib.finance import quotes_historical_yahoo
2+
from datetime import date
3+
import numpy as np
4+
import matplotlib.pyplot as plt
5+
from scipy import fftpack
6+
from scipy import signal
7+
from matplotlib.dates import DateFormatter, DayLocator, MonthLocator
8+
from scipy import optimize
9+
10+
11+
today = date.today()
12+
start = (today.year - 1, today.month, today.day)
13+
14+
quotes = quotes_historical_yahoo("QQQ", start, today)
15+
quotes = np.array(quotes)
16+
17+
dates = quotes.T[0]
18+
qqq = quotes.T[4]
19+
20+
21+
y = signal.detrend(qqq)
22+
23+
24+
alldays = DayLocator()
25+
months = MonthLocator()
26+
month_formatter = DateFormatter("%b %Y")
27+
28+
fig = plt.figure()
29+
ax = fig.add_subplot(211)
30+
31+
ax.xaxis.set_minor_locator(alldays)
32+
ax.xaxis.set_major_locator(months)
33+
ax.xaxis.set_major_formatter(month_formatter)
34+
35+
amps = np.abs(fftpack.fftshift(fftpack.rfft(y)))
36+
amps[amps < amps.max()] = 0
37+
38+
def residuals(p, y, x):
39+
A,k,theta = p
40+
err = y-A * np.sin(2* np.pi* k * x + theta)
41+
42+
return err
43+
44+
filtered = -fftpack.irfft(fftpack.ifftshift(amps))
45+
N = len(qqq)
46+
f = np.linspace(-N/2, N/2, N)
47+
p0 = [filtered.max(), f[amps.argmax()]/N, np.pi/3]
48+
plsq = optimize.leastsq(residuals, p0, args=(filtered, dates))
49+
p = plsq[0]
50+
print p
51+
plt.plot(dates, y, 'o', label="detrended")
52+
plt.plot(dates, filtered, label="filtered")
53+
plt.plot(dates, p[0] * np.sin(2 * np.pi * dates * p[1] + p[2]), '^', label="fit")
54+
fig.autofmt_xdate()
55+
plt.legend()
56+
57+
ax2 = fig.add_subplot(212)
58+
plt.plot(f, amps, label="transformed")
59+
60+
plt.legend()
61+
plt.show()

NumPy/Code/ch10code/frequencies.py

Lines changed: 54 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,54 @@
1+
from matplotlib.finance import quotes_historical_yahoo
2+
from datetime import date
3+
import numpy as np
4+
from scipy import signal
5+
import matplotlib.pyplot as plt
6+
from scipy import fftpack
7+
from matplotlib.dates import DateFormatter
8+
from matplotlib.dates import DayLocator
9+
from matplotlib.dates import MonthLocator
10+
11+
12+
today = date.today()
13+
start = (today.year - 1, today.month, today.day)
14+
15+
quotes = quotes_historical_yahoo("QQQ", start, today)
16+
quotes = np.array(quotes)
17+
18+
dates = quotes.T[0]
19+
qqq = quotes.T[4]
20+
21+
22+
y = signal.detrend(qqq)
23+
24+
25+
alldays = DayLocator()
26+
months = MonthLocator()
27+
month_formatter = DateFormatter("%b %Y")
28+
29+
fig = plt.figure()
30+
fig.subplots_adjust(hspace=.3)
31+
ax = fig.add_subplot(211)
32+
33+
ax.xaxis.set_minor_locator(alldays)
34+
ax.xaxis.set_major_locator(months)
35+
ax.xaxis.set_major_formatter(month_formatter)
36+
37+
# make font size bigger
38+
ax.tick_params(axis='both', which='major', labelsize='x-large')
39+
40+
amps = np.abs(fftpack.fftshift(fftpack.rfft(y)))
41+
amps[amps < 0.1 * amps.max()] = 0
42+
43+
plt.plot(dates, y, 'o', label="detrended")
44+
plt.plot(dates, -fftpack.irfft(fftpack.ifftshift(amps)), label="filtered")
45+
fig.autofmt_xdate()
46+
plt.legend(prop={'size':'x-large'})
47+
48+
ax2 = fig.add_subplot(212)
49+
ax2.tick_params(axis='both', which='major', labelsize='x-large')
50+
N = len(qqq)
51+
plt.plot(np.linspace(-N/2, N/2, N), amps, label="transformed")
52+
53+
plt.legend(prop={'size':'x-large'})
54+
plt.show()

NumPy/Code/ch10code/gaussquad.py

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
from scipy import integrate
2+
import numpy as np
3+
4+
print "Gaussian integral", np.sqrt(np.pi), integrate.quad(lambda x: np.exp(-x**2), -np.inf, np.inf)

NumPy/Code/ch10code/images.py

Lines changed: 31 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,31 @@
1+
from scipy import misc
2+
import numpy as np
3+
import matplotlib.pyplot as plt
4+
from scipy import ndimage
5+
6+
7+
image = misc.lena().astype(np.float32)
8+
9+
plt.subplot(221)
10+
plt.title("Original Image")
11+
img = plt.imshow(image, cmap=plt.cm.gray)
12+
plt.axis("off")
13+
14+
plt.subplot(222)
15+
plt.title("Median Filter")
16+
filtered = ndimage.median_filter(image, size=(42,42))
17+
plt.imshow(filtered, cmap=plt.cm.gray)
18+
plt.axis("off")
19+
20+
plt.subplot(223)
21+
plt.title("Rotated")
22+
rotated = ndimage.rotate(image, 90)
23+
plt.imshow(rotated, cmap=plt.cm.gray)
24+
plt.axis("off")
25+
26+
plt.subplot(224)
27+
plt.title("Prewitt Filter")
28+
filtered = ndimage.prewitt(image)
29+
plt.imshow(filtered, cmap=plt.cm.gray)
30+
plt.axis("off")
31+
plt.show()

NumPy/Code/ch10code/optfit.py

Lines changed: 67 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,67 @@
1+
from matplotlib.finance import quotes_historical_yahoo
2+
import numpy as np
3+
import matplotlib.pyplot as plt
4+
from scipy import fftpack
5+
from scipy import signal
6+
from matplotlib.dates import DateFormatter
7+
from matplotlib.dates import DayLocator
8+
from matplotlib.dates import MonthLocator
9+
from scipy import optimize
10+
11+
12+
start = (2010, 7, 25)
13+
end = (2011, 7, 25)
14+
15+
quotes = quotes_historical_yahoo("QQQ", start, end)
16+
quotes = np.array(quotes)
17+
18+
dates = quotes.T[0]
19+
qqq = quotes.T[4]
20+
21+
22+
y = signal.detrend(qqq)
23+
24+
25+
alldays = DayLocator()
26+
months = MonthLocator()
27+
month_formatter = DateFormatter("%b %Y")
28+
29+
fig = plt.figure()
30+
fig.subplots_adjust(hspace=.3)
31+
ax = fig.add_subplot(211)
32+
33+
ax.xaxis.set_minor_locator(alldays)
34+
ax.xaxis.set_major_locator(months)
35+
ax.xaxis.set_major_formatter(month_formatter)
36+
ax.tick_params(axis='both', which='major', labelsize='x-large')
37+
38+
amps = np.abs(fftpack.fftshift(fftpack.rfft(y)))
39+
amps[amps < amps.max()] = 0
40+
41+
def residuals(p, y, x):
42+
A,k,theta,b = p
43+
err = y-A * np.sin(2* np.pi* k * x + theta) + b
44+
45+
return err
46+
47+
filtered = -fftpack.irfft(fftpack.ifftshift(amps))
48+
N = len(qqq)
49+
f = np.linspace(-N/2, N/2, N)
50+
p0 = [filtered.max(), f[amps.argmax()]/(2*N), 0, 0]
51+
print "P0", p0
52+
53+
plsq = optimize.leastsq(residuals, p0, args=(filtered, dates))
54+
p = plsq[0]
55+
print "P", p
56+
plt.plot(dates, y, 'o', label="detrended")
57+
plt.plot(dates, filtered, label="filtered")
58+
plt.plot(dates, p[0] * np.sin(2 * np.pi * dates * p[1] + p[2]) + p[3], '^', label="fit")
59+
fig.autofmt_xdate()
60+
plt.legend(prop={'size':'x-large'})
61+
62+
ax2 = fig.add_subplot(212)
63+
ax2.tick_params(axis='both', which='major', labelsize='x-large')
64+
plt.plot(f, amps, label="transformed")
65+
66+
plt.legend(prop={'size':'x-large'})
67+
plt.show()

NumPy/Code/ch10code/pair.py

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,30 @@
1+
from matplotlib.finance import quotes_historical_yahoo
2+
from datetime import date
3+
import numpy as np
4+
from scipy import stats
5+
from statsmodels.stats.stattools import jarque_bera
6+
import matplotlib.pyplot as plt
7+
8+
9+
def get_close(symbol):
10+
today = date.today()
11+
start = (today.year - 1, today.month, today.day)
12+
13+
quotes = quotes_historical_yahoo(symbol, start, today)
14+
quotes = np.array(quotes)
15+
16+
return quotes.T[4]
17+
18+
spy = np.diff(np.log(get_close("SPY")))
19+
dia = np.diff(np.log(get_close("DIA")))
20+
21+
print "Means comparison", stats.ttest_ind(spy, dia)
22+
print "Kolmogorov smirnov test", stats.ks_2samp(spy, dia)
23+
24+
print "Jarque Bera test", jarque_bera(spy - dia)[1]
25+
26+
plt.hist(spy, histtype="step", lw=1, label="SPY")
27+
plt.hist(dia, histtype="step", lw=2, label="DIA")
28+
plt.hist(spy - dia, histtype="step", lw=3, label="Delta")
29+
plt.legend()
30+
plt.show()
Lines changed: 31 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,31 @@
1+
from scipy.io import wavfile
2+
import matplotlib.pyplot as plt
3+
import urllib2
4+
import numpy as np
5+
import sys
6+
7+
response = urllib2.urlopen('http://www.thesoundarchive.com/austinpowers/smashingbaby.wav')
8+
print response.info()
9+
WAV_FILE = 'smashingbaby.wav'
10+
filehandle = open(WAV_FILE, 'w')
11+
filehandle.write(response.read())
12+
filehandle.close()
13+
sample_rate, data = wavfile.read(WAV_FILE)
14+
print "Data type", data.dtype, "Shape", data.shape
15+
16+
plt.subplot(2, 1, 1)
17+
plt.title("Original")
18+
plt.plot(data)
19+
20+
plt.subplot(2, 1, 2)
21+
22+
# Repeat the audio fragment
23+
repeated = np.tile(data, int(sys.argv[1]))
24+
25+
# Plot the audio data
26+
plt.title("Repeated")
27+
plt.plot(repeated)
28+
wavfile.write("repeated_yababy.wav",
29+
sample_rate, repeated)
30+
31+
plt.show()

NumPy/Code/ch10code/scipyio.py

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,6 @@
1+
import numpy as np
2+
from scipy import io
3+
4+
a = np.arange(7)
5+
6+
io.savemat("a.mat", {"array": a})

0 commit comments

Comments
 (0)