2013-01-07 3 views
2

나는 다음 배포 한 :R의 곡선 아래 면적의 95 % 신뢰 한계를 계산하려면 어떻게해야합니까?

x<-c(22.5,28.14285714,33.78571429,39.42857143,45.07142857,50.71428571,56.35714286,62,67.64285714,73.28571429,78.92857143,84.57142857,90.21428571,95.85714286,101.5,107.1428571,112.7857143,118.4285714,124.0714286,129.7142857,135.3571429,141,146.6428571,152.2857143,157.9285714,163.5714286,169.2142857,174.8571429,180.5,186.1428571,191.7857143,197.4285714,203.0714286,208.7142857,214.3571429,220,225.6428571,231.2857143,236.9285714,242.5714286,248.2142857,253.8571429,259.5,265.1428571,270.7857143,276.4285714,282.0714286,287.7142857,293.3571429,299) 
y<-c(0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.00328839614285714,0.00296425985714286,0.002655899,0.00236187857142857,0.002080895,0.00181184271428571,0.00155376085714286,0.00130578928571429,0.001074706,0.000877193,0.000709397142857142,0.000567189714285714,0.000447254,0.000346858571428571,0.000263689142857143,0.000195768428571429,0.000141427,9.92657142857141e-05,6.77857142857142e-05,4.48571428571428e-05,2.86428571428571e-05,1.75142857142857e-05,1.01357142857143e-05,5.52e-06,2.78857142857142e-06,1.27285714285713e-06,5.00714285714284e-07,1.5742857142857e-07,3.29857142857142e-08,2.78857142857137e-09,1.74e-12) 

plot(x,y) 

내가 왼쪽과 오른쪽 영역의 0.05 분포에서 0.95의 영역을 분리 x의 가치를 발견하고 싶습니다 (95 %를 하나의 꼬리 신뢰성의 간격).

필자는 경험적 곡선을 함수에 적용한 다음 원하는 값을 얻을 수 있도록 함수를 통합해야하지만 실제로 어디서부터 시작해야할지 모릅니다.

어떻게 R에서이 작업을 수행 할 수 있습니까?

+0

GSee의 대답은 갈 길입니다. 그러나 원본 데이터의 수치 적 통합은 적합 함수를 작성하고 통합하는 것보다 쉬울뿐만 아니라 일반적으로 계산 오류가 적음을 지적하고자합니다. –

+2

@CarlWitthoft, 나는 (quantile (x, 0.95)) 내 대답에 대해서는 그렇게 확신하지 못한다. 'x'를 95 %와 5 %로 나누지 만, 영역 (y's)는 전혀 고려하지 않습니다. – GSee

+0

@Gsee - 결국 심슨의 적분 값을 생성해야한다고 생각합니다. 나는 여전히 fit 함수를 생성하지 않는다. –

답변

2

통합 문제 (곡선 아래의 합계)입니다. 통합을 정사각형 + 곡선으로 나눌 수 있습니다. 그러나, 당신은 스플라인을 통해 신속하고 더러운 근사치를 사용할 수 있습니다

y<-c(0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.003541755,0.00328839614285714,0.00296425985714286,0.002655899,0.00236187857142857,0.002080895,0.00181184271428571,0.00155376085714286,0.00130578928571429,0.001074706,0.000877193,0.000709397142857142,0.000567189714285714,0.000447254,0.000346858571428571,0.000263689142857143,0.000195768428571429,0.000141427,9.92657142857141e-05,6.77857142857142e-05,4.48571428571428e-05,2.86428571428571e-05,1.75142857142857e-05,1.01357142857143e-05,5.52e-06,2.78857142857142e-06,1.27285714285713e-06,5.00714285714284e-07,1.5742857142857e-07,3.29857142857142e-08,2.78857142857137e-09,1.74e-12) 
x<-c(22.5,28.14285714,33.78571429,39.42857143,45.07142857,50.71428571,56.35714286,62,67.64285714,73.28571429,78.92857143,84.57142857,90.21428571,95.85714286,101.5,107.1428571,112.7857143,118.4285714,124.0714286,129.7142857,135.3571429,141,146.6428571,152.2857143,157.9285714,163.5714286,169.2142857,174.8571429,180.5,186.1428571,191.7857143,197.4285714,203.0714286,208.7142857,214.3571429,220,225.6428571,231.2857143,236.9285714,242.5714286,248.2142857,253.8571429,259.5,265.1428571,270.7857143,276.4285714,282.0714286,287.7142857,293.3571429,299) 

sp=smooth.spline(x,y) 
f = function(t) 
{ 
    predict(sp,t)$y 
} 

N=500 # this is an accuracy parameter 
xBis=seq(x[1],x[length(x)],length=N) 
yBis=sapply(x,f) 

J = function (input) 
{ # This function takes input in 1:N 
    Integral = 0 
    dx=(x[length(x)]-x[1])/N 

    for (j in 1: input) 
{ z=xBis[j] 
    Integral=Integral+ f(x[1]+z)*dx 
} 
J=Integral 
} 
###### 
I=J(N) # This is the value of the sum under the curve 
# It should be roughly equal (given the shape of the curve) to: 
index=max(which(y==y[1])) 
I = (x[index]-x[1])*(y[index])*3/2 
###### 
res=sapply(1:N,J)/I 
Index5=max(which(res<=.05)) 
Index95=min(which(res>=.95)) 

x5=xBis[Index5] # This is the 5% quantile 
x95=xBis[Index95] 

HTH

아무것도 내가 할 수있는 훨씬 더 나은 방법이 생각

PS 불분명 한 경우 알려주세요 ..

4

다른 답변에서 지적한 바와 같이, 이것은 곡선 문제에서 통합되어 전체 면적의 95 %에 도달하는 위치와 결부됩니다. 나는 David's answer보다 단순한 통합 접근법을 택했다. 커브를 보간하고 통합하는 대신 사다리꼴 통합 규칙을 사용하여 각 간격으로 제공된 영역을 얻습니다. 그런 다음 개별 영역을 추가하여 전체 영역을 가져옵니다. 누적 영역이 전체 영역의 95 %를 초과하는 색인이 발견되고이를 사용하여 회선을 그릴 수 있습니다.

piece_area <- c(0, (x[-1] - x[-length(x)])*(y[-1] + y[-length(y)])/2) 
cum_area <- cumsum(piece_area) 
total_area <- cum_area[length(cum_area)] 
idx095 <- min(which(cum_area > 0.95 * total_area)) 

abline(v = x[idx095]) 
95 %가 분포의 원래 샘플에서 더 많은 포인트를 사용하여 얻을 수 교차되는 정확한 지점의

enter image description here

높은 해상도.

관련 문제