《天文数据分析课件.ppt》由会员分享，可在线阅读，更多相关《天文数据分析课件.ppt（50页珍藏版）》请在三一办公上搜索。
1、天文数据分析,国家天文台 赵永恒2015年4月,数据分析方法,累积概率分布,任意分布的随机数,常用概率分布,数据插值,内插和外推,插值(scipy.interpolate),线性1d插值(interp1d)1d样条插值(interpolate.splXXX)2d样条插值(bisplrep)RBF（radial basis function）插值/平滑,import numpy as npfrom scipy.interpolate import Rbfimport matplotlib.pyplot as pltfrom matplotlib import cm#2-d tests-setup
2、 scattered datax=np.random.rand(100)*4.0-2.0y=np.random.rand(100)*4.0-2.0z=x*np.exp(-x*2-y*2)ti=np.linspace(-2.0,2.0,100)XI,YI=np.meshgrid(ti,ti)#use RBFrbf=Rbf(x,y,z,epsilon=2)ZI=rbf(XI,YI)#plot the resultn=plt.Normalize(-2.,2.)plt.subplot(1,1,1)plt.pcolor(XI,YI,ZI,cmap=cm.jet)plt.scatter(x,y,100,z
3、,cmap=cm.jet)plt.title(RBF interpolation-multiquadrics)plt.xlim(-2,2)plt.ylim(-2,2)plt.colorbar()plt.show(),数据拟合,拟合模型,O=T+e,最小二乘法,import numpy as npimport scipy.linalg as splinimport matplotlib.pyplot as pltc1,c2=5.0,2.0i=np.r_1:11xi=0.1*iyi=c1*np.exp(-xi)+c2*xizi=yi+0.05*np.max(yi)*np.random.randn(
4、len(yi)A=np.c_np.exp(-xi):,np.newaxis,xi:,np.newaxisc,resid,rank,sigma=splin.lstsq(A,zi)print c,resid,rank,sigmaxi2=np.r_0.1:1.0:100jyi2=c0*np.exp(-xi2)+c1*xi2plt.plot(xi,zi,x,xi2,yi2)plt.axis(0,1.1,3.0,5.5)plt.xlabel($x_i$)plt.title(Data fitting with linalg.lstsq)plt.show(),非线性拟合,拟合系数有非线性函数关系,from
5、numpy import x=arange(0,6e-2,6e-2/30)A,k,theta=10,1.0/3e-2,pi/6y_true=A*sin(2*pi*k*x+theta)y_meas=y_true+2*random.randn(len(x)def residuals(p,y,x):A,k,theta=p err=y-A*sin(2*pi*k*x+theta)return errdef peval(x,p):return p0*sin(2*pi*p1*x+p2)p0=8,1/2.3e-2,pi/3print array(p0)#8.43.4783 1.0472from scipy.o
6、ptimize import leastsqplsq=leastsq(residuals,p0,args=(y_meas,x)print plsq0#10.9437 33.3605 0.5834print array(A,k,theta)#10.33.3333 0.5236import matplotlib.pyplot as pltplt.plot(x,peval(x,plsq0),x,y_meas,o,x,y_true)plt.title(Least-squares fit to noisy data)plt.legend(Fit,Noisy,True)plt.show(),Maximum
7、 Likelihood Estimation(MLE),最大似然估计,一致收敛渐近于正态分布：有极值,等精度高斯分布,非等精度高斯分布,上下限分布,求极值scipy.optimize,A collection of general-purpose optimization routines.fmin-Nelder-Mead Simplex algorithm(uses only function calls)fmin_powell-Powells(modified)level set method(uses only function calls)fmin_cg-Non-linear(Pola
8、k-Ribiere)conjugate gradient algorithm(can use function and gradient).fmin_bfgs-Quasi-Newton method(Broydon-Fletcher-Goldfarb-Shanno);(can use function and gradient)fmin_ncg-Line-search Newton Conjugate Gradient(can use function,gradient and Hessian).leastsq-Minimize the sum of squares of M equation
9、s in N unknowns given a starting estimate.Constrained Optimizers(multivariate)fmin_l_bfgs_b-Zhu,Byrd,and Nocedals L-BFGS-B constrained optimizer(if you use this please quote their papers-see help)fmin_tnc-Truncated Newton Code originally written by Stephen Nash and adapted to C by Jean-Sebastien Roy
10、.fmin_cobyla-Constrained Optimization BY Linear Approximation,EM算法,Expect-Maximum算法EM算法的目标是找出有隐性变量的概率模型的最大可能性解分为2个步骤，E-step和M-stepE-step根据最初假设的模型参数值或者上一步的模型参数计算出隐性变量的后验概率，其实就是隐性变量的期望M-step根据这个E-step的后验概率重新计算出模型参数，然后再重复这两个步骤，直至目标函数收敛。,拟合度,好模型,精度估计,Bootstrap方法“杰克刀”方法,Bootstrap方法,用已知数据建立分布可有N!组新数据,SN19
11、87A中微子事件,9个中微子,“杰克刀”方法,每次丢掉1或多个数据,KolmogorovSmirnov(K-S)检验,高斯分布？,AndersonDarling检验ShapiroWilk检验,Bayes统计推断,高斯分布,问题：,问题：,直方图,Markov chain Monte Carlo(MCMC),https:/,时间序列分析,Fourier变换scipy.fftpack,卷积定理,离散Fourier变换,窗效应,GX 5-1,LombScargle 周期图,不等间隔,相位图,注意：倍频、差频,相关分析,小波变换,课后作业II,1、一组有确定均值和误差的测量值中，即有N个值xi，又有M个上下限值xmini,xmaxi，请写出其似然函数ln L。2、随机写出20个整数（可重复），组成xi。以xi的分布产生出1000个随机数，构成yi。并用KS检验说明xi和yi的分布是否一致。3、设计具有一个或数个频率（周期）的正弦曲线，并加上误差。用Fourier谱分析确定其频率（周期）。,
