1. 程式人生 > >【高階程式設計技術】第12周作業

【高階程式設計技術】第12周作業

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as sts

#exercise 1

f, ax = plt.subplots(1, 1, figsize=(5,4))

x=np.linspace(0,2,50)
y=np.multiply(np.power(np.sin(np.subtract(x,2)),2),np.exp(np.negative(np.power(x,2))))

ax.plot(x, y)
ax.set_xlim((0,2))

ax.set_xlabel('my x label')
ax.set_ylabel('my y label')
ax.set_title('plot title')


plt.show()


import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import leastsq

data = np.array([])
for sigma in range(1,11):
	data=np.append(data,sigma*np.random.randn(1,20) + sigma*sigma)

data=np.reshape(data,(10,20))
data=data.T
#20*10

b=np.random.rand(10,1)#10*1
z=np.random.randn(20,1)#20*1
y=np.add(np.dot(data,b),z)#20*1

def fun(p,X):
	b=np.reshape(p,(10,1))
	return np.dot(X,b).flatten()


def err(p,X,y):
	return np.subtract(fun(p,X),y)

p0=np.random.rand(10,1)

y=y.flatten()#將y變為一維陣列
res=leastsq(err,p0,args=(data,y))

b=b.flatten()#將b變為一維陣列

x_ax=np.arange(0,10)#橫軸為0到10
y_ax=res[0][x_ax]#估計出的b
z_ax=b[x_ax]#隨機生成的b
f, ax = plt.subplots(1, 1, figsize=(5,4))
ax.scatter(x_ax,y_ax,c='r',marker='o')
ax.scatter(x_ax,z_ax,c='b',marker='x')
plt.show()


import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as sts

#exercise 3

data = np.random.normal(size=(10000,1))
data_mean=np.mean(data,axis=1)

kde=sts.gaussian_kde(data_mean)

x=np.linspace(-10,10,50)

f, ax= plt.subplots(1,2,figsize=(6,3))
ax[0].hist(data,bins=25,color='b')
ax[1].plot(x,kde(x),label='kde',color='r')

#print(res)
plt.show()