我需要将一些理论数据与
python中的实际数据进行比较.
理论数据来自解决方程式.
为了改进比较,我想删除远离理论曲线的数据点.我的意思是,我想删除图中红色虚线下方和上方的点(用matplotlib制作).
理论数据来自解决方程式.
为了改进比较,我想删除远离理论曲线的数据点.我的意思是,我想删除图中红色虚线下方和上方的点(用matplotlib制作).
理论曲线和数据点都是不同长度的阵列.
我可以尝试以粗略的方式移除点,例如:可以使用以下方法检测第一个上点:
data2[(data2.redshift<0.4)&data2.dmodulus>1] rec.array([('1997o',0.374,1.0203223485103787,0.44354759972859786)],dtype=[('SN_name','|S10'),('redshift','<f8'),('dmodulus',('dmodulus_error','<f8')])
但我想用一种不太粗略的方式.
那么,任何人都可以帮我找到一个简单的方法来消除问题点吗?
谢谢!
解决方法
这可能是过度的,并且基于您的评论
Both the theoretical curves and the data points are arrays of
different length.
我会做以下事情:
>截断数据集,使其x值位于理论集的最大值和最小值之内.
>使用scipy.interpolate.interp1d和上面截断的数据x值插值理论曲线.步骤(1)的原因是满足interp1d的约束.
>使用numpy.where查找超出可接受理论值范围的数据y值.
>不要像评论和其他答案中所建议的那样丢弃这些价值观.如果你想要清晰,可以通过绘制’内衬’一种颜色和’异常值’另一种颜色来指出它们.
我认为这是一个接近你所寻找的脚本.它有望帮助您实现您想要的目标:
import numpy as np import scipy.interpolate as interpolate import matplotlib.pyplot as plt # make up data def makeUpData(): '''Make many more data points (x,y,yerr) than theory (x,y),with theory yerr corresponding to a constant "sigma" in y,about x,y value''' NX= 150 datax = (np.random.rand(NX)*1.1)**2 dataY = (1.5*datax+np.random.rand(NX)**2)*datax dataErr = np.random.rand(NX)*datax*1.3 theoryX = np.arange(0,1,0.1) theoryY = theoryX*theoryX*1.5 theoryErr = 0.5 return datax,dataY,dataErr,theoryX,theoryY,theoryErr def makeSameXrange(theoryX,datax,dataY): ''' Truncate the datax and dataY ranges so that datax min and max are with in the max and min of theoryX. ''' minT,maxT = theoryX.min(),theoryX.max() goodIdxMax = np.where(datax<maxT) goodIdxMin = np.where(datax[goodIdxMax]>minT) return (datax[goodIdxMax])[goodIdxMin],(dataY[goodIdxMax])[goodIdxMin] # take 'theory' and get values at every 'data' x point def theoryYatdatax(theoryX,datax): '''For every datax point,find interpolated thoeryY value. theoryx needed for interpolation.''' f = interpolate.interp1d(theoryX,theoryY) return f(datax[np.where(datax<np.max(theoryX))]) # collect valid points def findInlierSet(datax,interpTheoryY,thoeryErr): '''Find where theoryY-theoryErr < dataY theoryY+theoryErr and return valid indicies.''' withinUpper = np.where(dataY<(interpTheoryY+theoryErr)) withinLower = np.where(dataY[withinUpper] >(interpTheoryY[withinUpper]-theoryErr)) return (datax[withinUpper])[withinLower],(dataY[withinUpper])[withinLower] def findOutlierSet(datax,thoeryErr): '''Find where theoryY-theoryErr < dataY theoryY+theoryErr and return valid indicies.''' withinUpper = np.where(dataY>(interpTheoryY+theoryErr)) withinLower = np.where(dataY<(interpTheoryY-theoryErr)) return (datax[withinUpper],dataY[withinUpper],datax[withinLower],dataY[withinLower]) if __name__ == "__main__": datax,theoryErr = makeUpData() Truncdatax,TruncDataY = makeSameXrange(theoryX,dataY) interpTheoryY = theoryYatdatax(theoryX,Truncdatax) indatax,inDataY = findInlierSet(Truncdatax,TruncDataY,theoryErr) outUpX,outUpY,outDownX,outDownY = findOutlierSet(Truncdatax,theoryErr) #print inlierIndex fig = plt.figure() ax = fig.add_subplot(211) ax.errorbar(datax,fmt='.',color='k') ax.plot(theoryX,'r-') ax.plot(theoryX,theoryY+theoryErr,'r--') ax.plot(theoryX,theoryY-theoryErr,'r--') ax.set_xlim(0,1.4) ax.set_ylim(-.5,3) ax = fig.add_subplot(212) ax.plot(indatax,inDataY,'ko') ax.plot(outUpX,'bo') ax.plot(outDownX,outDownY,'ro') ax.plot(theoryX,3) fig.savefig('findInliers.png')
这个数字是结果: