2013-06-07 10 views
3

데이터 프레임의 모든 고유 사이트에 대한 풀 (실제로 가중치 적용) 표준 편차를 계산하고 싶습니다.R에서 풀 표준 편차를 계산하는 방법?

이 사이트의 값은 단일 종 숲 스탠드의 값이며 평균값과 sd를 모아서 침착 물 스탠드와 폭 넓은 스탠드를 비교할 수 있습니다.
이것은 된 Broadleaved 스탠드 값과 데이터 프레임 (DF)가있다 :

keybl   n mean sd 
Vest02DenmDesp 3 58.16 6.16 
Vest02DenmDesp 5 54.45 7.85 
Vest02DenmDesp 3 51.34 1.71 
Vest02DenmDesp 3 59.57 5.11 
Vest02DenmDesp 5 62.89 10.26 
Vest02DenmDesp 3 77.33 2.14 
Mato10GermDesp 4 41.89 12.6 
Mato10GermDesp 4 11.92 1.8 
Wawa07ChinDesp 18 0.097 0.004 
Chen12ChinDesp 3 41.93 1.12 
Hans11SwedDesp 2 1406.2 679.46 
Hans11SwedDesp 2 1156.2 464.07 
Hans11SwedDesp 2 4945.3 364.58 

Keybl 사이트에 대한 코드이다. 풀링 된 SD에 대한 공식은 다음과 같습니다

s=sqrt((n1-1)*s1^2+(n2-1)*s2^2)/(n1+n2-2)) 

2

그룹의 수이고, 따라서 (내가 사진을 게시 할 수 없습니다 직접 공식에 갈 것 링크를 찾을 수 없습니다 죄송합니다) 사이트에 따라 변경됩니다. 나는 이것이 t-test에 사용되고 두 그룹이 비교하기를 원한다는 것을 안다. 이 경우 나는 이들 그룹을 비교할 계획이 아닙니다. 교수님은이 공식을 사용하여 가중치 적용 SD를 얻으라고 제안했습니다. 나는 내가 필요로하는 방식으로이 공식을 통합 한 R 함수를 찾지 못했고, 그래서 나는 내 자신을 만들려고 노력했다. 그러나 저는 R에 익숙하지 않고 함수와 루프를 만드는 데별로 능숙하지 않으므로 여러분의 도움을 바랍니다.

이것은 내가 지금까지 무엇을 가지고 있습니다 : 그것은 약간 낮은 값이해야 내가 왜 이해하지 못하는 이상을 제공으로

sd=function (data) { 
nc1=data[z,"nc"] 
sc1=data[z, "sc"] 
nc2=data[z+1, "nc"] 
sc2=data[z+1, "sc"] 
sd1=(nc1-1)*sc1^2 + (nc2-1)*sc2^2 
sd2=sd1/(nc1+nc2-length(nc1)) 
sqrt(sd2) 
} 

splitdf=split(df, with(df, df$keybl), drop = TRUE) 

for (c in 1:length(splitdf)) { 
for (i in 1:length(splitdf[[i]])) { 
    a = (splitdf[[i]]) 
    b =sd(a) 
    } 
} 

1) 함수 자체가 올바르지 않습니다. z + 1이 마지막 행에 도달했을 때 멈추지 않을 수 있습니까? 그렇다면 어떻게 수정 될 수 있습니까?

2) 루프가 완전히 잘못되었지만 몇 시간 동안 성공하지 못했을 때 다시 생각해 낼 수 있습니다.

아무도 도와 줄 수 있습니까?

감사합니다,

Antra 당신이 그것을 쉽게 할 것이다보다 일반적인 식에 도움이 될 뭘 하려는지

+0

: S = SQRT의 (((N1-1) S1 *^2 + (n2-1) S2 *^2)/(N1 + n2-2)) –

답변

2

독립성 가정하에 (따라서 공분산 항은 0으로 가정 할 수 있습니다) 풀링 된 SD는 다음과 같습니다. sqrt (sum_over_groups [(var)/sum (n) -N_groups)])

여기서 괄호 누락
 lapply(split(dat, dat$keybl), 
      function(dd) sqrt(sum(dd$sd^2 * (dd$n-1))/(sum(dd$n-1)-nrow(dd)))) 
#------------------------- 
$Chen12ChinDesp 
[1] 1.583919 

$Hans11SwedDesp 
[1] Inf 

$Mato10GermDesp 
[1] 11.0227 

$Vest02DenmDesp 
[1] 9.003795 

$Wawa07ChinDesp 
[1] 0.004123106 
+0

'lapply + split' ~'by'? ;-) – agstudy

+0

상호 이변 조건을 0으로 가정해야한다는 점을 인정하지 않습니까? 요점은 : 때로는'do.call (rbind (.)) '이 함께 되돌아 가야하고 clunky (두 개가 아닌 세 개의 함수)가 필요하기 때문에 피할 수 있지만'sapply (split)) '여기에 있습니다. 내가 생각하는 스타일의 문제. –

+1

고마워요. 이 대답은 잘 작동하고 내가 생각했던 것과 비교할 때 훨씬 좋습니다. 나는 sd가 n-1로 곱 해져야하므로 (dd $ n-1) 함수에 추가했다. sum (dd $ sd^2 * (dd $ n-1)) – Antra

5

. keybl 변수에 의해 조각으로 나눌 필요가 없다면 할 수 있습니다.

dd <- df #df is not a good name for a data.frame variable since df has a meaning in statistics 

dd$df <- dd$n-1 
pooledSD <- sqrt(sum(dd$sd^2 * dd$df)/sum(dd$df)) 
# note, in this case I only pre-calculated df because I'll need it more than once. The sum of squares, variance, etc. are only used once. 

R의 중요한 일반 원칙은 가능한 한 벡터 수학을 사용한다는 것입니다. 이 사소한 경우에는별로 중요하지 않지만 계산 속도가 더 중요하다고 생각되는 큰 data.frame 개체에서 이것을 수행하는 방법을 확인하려면 읽으십시오.

# First use R's vector facilities to define the variables you need for pooling. 
dd$df <- dd$n-1 
dd$s2 <- dd$sd^2 # sd isn't a good name for standard deviation variable even in a data.frame just because it's a bad habit to have... it's already a function and standard deviations have a standard name 
dd$ss <- dd$s2 * dd$df 

이제는 필요한 합계를 나누고 계산하기위한 편리한 함수를 사용하십시오. 암시 적 루프마다 하나의 함수 만 실행됩니다 (* apply, aggregate 등은 함수를 여러 번 실행하는 암시 적 루프입니다).

ds <- aggregate(ss ~ keybl, data = dd, sum) 
ds$df <- tapply(dd$df, dd$keybl, sum) #two different built in methods for split apply, we could use aggregate for both if we wanted 
# divide your ss by your df and voila 
ds$s2 <- ds$ss/ds$df 
# and also you can easly get your sd 
ds$s <- sqrt(ds$s2) 

그리고 정답은 :

  keybl   ss df   s2   s 
1 Chen12ChinDesp 2.508800e+00 2 1.254400e+00 1.120000 
2 Hans11SwedDesp 8.099454e+05 3 2.699818e+05 519.597740 
3 Mato10GermDesp 4.860000e+02 6 8.100000e+01 9.000000 
4 Vest02DenmDesp 8.106832e+02 16 5.066770e+01 7.118125 
5 Wawa07ChinDesp 2.720000e-04 17 1.600000e-05 0.004000 

이 훨씬 덜 간결 다른 방법보다 본다 (42의 대답처럼)하지만 당신은 얼마나 많은 R 명령의 관점에서 사람들을 풀다 경우 실제로 이것을 실행하는 것이 훨씬 간결합니다.이와 같은 짧은 문제에 대해서는 어느 쪽이든 괜찮 으면하지만 대부분의 벡터 수학을 사용하는 방법을 보여줄 것이라고 생각했습니다. 또한 표현의 편의를 위해 이러한 편리한 암시 적 루프 기능을 사용할 수있는 이유를 강조합니다. 만약 당신이 for 루프를 사용했다면 같은 것을 이루려면 모든 것을 루프에 넣으려고하는 것이 더 강합니다. R에서 나쁜 생각 일 수 있습니다.

+1

비행 중에'n-1' 만 하는게 어떨까요? 'sqrt (sum (df $ sd^2 * (df $ n - 1))/(합계 (df $ n - 1)))) ' – atomicules

+0

좀 더 명확하게 대답을 확장했습니다. – John

관련 문제