#acl +All:read
#format wiki
#language ko
#pragma description 기초의학통계학 및 실험;
= R Basic =

== 기본 of 기본 ==
R을 설치하면 가지고 놀 수 있는 [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/00Index.html|datasets]]가 기본으로 설치된다. [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/00Index.html|한번 들여다보자]]
 * [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/trees.html|trees]] (벗나무의 직경, 높이, 부피)
 * [[WikiPedia:Iris flower data set]] (붓꽃 길이)
 * [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/ChickWeight.html|ChickWeight]] (병아리 몸무게)
 * [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/Loblolly.html|Loblolly]] (Loblolly 소나무 높이)
 * [[https://stat.ethz.ch/R-manual/R-devel/library/datasets/html/Orange.html|Orange]] (오렌지 나무 직경)

=== trees dataset ===
 * 벗나무 31그루의 직경, 높이, 부피를 구한 데이타이다.
{{{#!highlight r
# trees란 무엇인가!
?trees

# trees 데이타를 화면에 띄워보자
trees # 혹은 print(trees)
str(trees)
table(trees)

# data size(dimension)을 확인해보자.
dim(trees)

# 31 x 3 개의 숫자를 한번에 보려니 눈에도 안 들어오고 의미없다. 일부만 화면에 띄워보자.
head(trees) # 첫 6개의 data point만 화면에 출력
head(trees,10) # 첫 10개의 data point만 화면에 출력
tail(trees,4) # 마지막 4개의 data point만 화면에 출력

# 데이타의 평균, 최상값, 최소값, 중간값, 상위 25% 값, 하위 25% 값을 보고 싶다면,
summary(trees)

# 나무 높이의 평균, 최상값, 최소값, 중간값, 상위 25% 값, 하위 25% 값을 따로 보고 싶다면,
mean(trees$Height) # mean(trees) 은 어떤 결과가 나오나?
max(trees$Height) # max(trees) 은 어떤 결과가 나오나?
min(trees$Height) # min(trees) 은 어떤 결과가 나오나?
quantile(trees$Height) # quantile(trees) 은 어떤 결과가 나오나?
}}}

숫자보다 그림이 데이타를 한눈에 파악하기 편하다.
{{{#!highlight r
boxplot(trees)

# 통계 계산값 말고, 모든 data point를 다 xy 평면에 찍어줘
plot(trees)

# 특별히 Girth 대비 Volume만 따로 xy 평면에 보여주되, xy 둘다 log scale 로 보겠어
plot(Volume ~ Girth, data = trees, log = "xy")
# trees 데이타의 컬럼 이름을 지정하는 방법
plot(trees$Girth, trees$Volume, log = "xy")
# trees 데이타의 컬럼 번호를 지정하는 방법
plot(trees[,1], trees[,3], log = "xy")

trees[c(3:10), c(1,3)]
plot(trees[c(3:10), c(1,3)])
}}}
Tip:
 * data[row 번호, col 번호] 와 같이 data point의 index로 지정가능
 * 여러개의 숫자(혹은 문자)를 묶을 때는 c( , , ) 혹은 c( : ) 형식을 사용함

궁극의 Tip:
 * Console에서 위화살표 키를 누르면 전에 쳤던 명령문이 복사된다.
 * #이후의 문장은 실행이 안 된다. 메모/도움말 용도로 사용하자

=== iris dataset ===
 * [[WikiPedia:Iris (plant)|붓꽃]] 150송이의 꽃받침 길이(Sepal.Length), 꽃받침 너비(Sepal.Width), 꽃잎 길이(Petal.Length), 꽃잎 너비(Petal.Width), 종(Species)을 조사한 데이터이다.
 * 생물학자&통계학자인 [[WikiPedia:Ronald Fisher]]의 1936년 [[WikiPedia:Linear discriminant analysis]]에 관한 논문에 등장한 이후로 [[WikiPedia:Multivariate statistics]]의 표준 예제 데이타로 자주 인용된다.

{{{#!highlight r
# 다음을 실행 해 보자
iris
dim(iris) # data size(dimension)을 확인해보자.
head(iris) # 첫 6개의 data point만 화면에 출력
tail(iris,5) # 마지막 5개의 data point만 화면에 출력
summary(iris)
boxplot(iris)
plot(iris) # scatter plot

# Species에 따라 다른 색을 칠하면 보기 편하다
plot(iris, col=c('red','blue','green3')[unclass(iris$Species)])

# 특정 측정값, 예를 들어 Petal.Length vs. Sepal.Length 만 보고 싶다면,
plot(Sepal.Length ~ Petal.Length, data = iris, col=c('red','blue','green3')[unclass(iris$Species)] )
plot(iris$Petal.Length, iris$Sepal.Length)

# 한동안 iris 데이타 가지고 놀 것 같으면 미리 선언할 수 있다.
attach(iris) # 앞으로 언급 안 하면 data = iris 라는 뜻임
plot(Petal.Length, Sepal.Length, col=c('red','blue','green4')[unclass(Species)])
plot(Sepal.Length ~ Petal.Length, col=c('red','blue','green4')[unclass(Species)])
detach(iris) # iris 놀이 끝
}}}

=== ChickWeight ===
병아리 50마리를 대상으로 조사한 다이어트별 성장 속도를 조사한 데이타이다.
 * weight: 병아리 몸무게 (numeric 데이타, gram)
 * Time: 나이 (numeric 데이타, days)
 * Chick: 병아리 일련번호 (ordered factor)
 * Diet: 다이어트 종류 (factor)

{{{
과제: ChickWeight 데이타를 이용하여, 다이어트 종류별로 다른 색을 사용하여, 시간 vs 몸무게 그래프를 그려라.
}}}

<<Include(../MidTerm)>>

== Birthday Paradox ==
[[WikiPedia:Birthday Paradox]]: 23명이 모이면 그 안에 생일이 같은 사람이 존재할 확률은 50% 이상이다.
 * 2014 FIFA 월드컵 본선 진출 32팀(각 팀 23명) 중 생일이 같은 멤버가 있던 팀은 16팀이었음 ([[http://www.bbc.co.uk/news/magazine-27835311|BBC]])
 * [[http://betterexplained.com/articles/understanding-the-birthday-paradox/|참고자료]]

계산해보자
 * 365개의 공을 중복 가능하게 뽑을 때, 같은 공이 나올 확률 문제이다.
 * 같은 공이 나올 확률은 (1 - 다른 공이 나올 확률) 이다.
 * 다른 공이 나올 확률은 횟수를 거듭할수록 1/365 씩 줄어든다.
 * n 번째까지 모두 다른 공이 나올 확률은 각 회당 다른 공이 나올 확률의 누적곱이다.
{{{#!highlight r
person <- 23
odd <- 1.0

for (n in 2:person) {
  odd = odd * (366 - n) / 365
  print((1-odd) * 100)
}
}}}

참고로 python이라면 이렇게 쓴다
{{{#!highlight python
person = 23
odd = 1.0

for n in range(2,person+1):
    odd = odd * (366 - n)/365
    print (1-odd) * 100
}}}

수학이 어려우면 무식한 시뮬레이션으로 알아볼 수도 있다 (brute force algorithm)
 * 함수를 하나 만들어서 birthday라고 이름을 붙힌다.
 * birthday 함수는 nPeople, num 두 가지 변수를 input으로 받는다.
 * input value를 생략하면 두 변수의 default값은 각각 23, 1000이다.
 * 함수가 하는 일은 1에서 365 사이에서 중복허용으로 nPeople 개의 숫자를 고른 후, 그 중 중복이 있는지를 판단하는 것을 num번 반복하여, 그 중 몇 %에서 중복이 있는지 출력한다.
{{{#!highlight r
birthday <- function(nPeople = 23, num = 1000) {
  data <- rep(NA,num)
  for (i in 1:num) {
    thisGroup <- sample(1:365, nPeople, replace = TRUE)
    data[i] <- sum(duplicated(thisGroup)) > 0
  }
  print(sum(data)/num)
}
birthday()
birthday(nPeople = 30, num =10000)
}}}

for-loop으로 구성된 코드는 가독성(human readability)은 좋지만 처리속도는 불리하다.
 * R에서는 for-loop대신 sapply 함수가 선호된다. sapply(이 숫자들을, function(x)에 넣어라) 라는 뜻이다.
{{{#!highlight r
birthdays <- sapply(1:10000, function(x){sample(1:365, 23, replace = TRUE)})
sum(sapply(1:10000, function(x) sum(duplicated(birthdays[,x])) > 0))

# nice one-liner
birthdays <- function(x){sum(duplicated(sample(1:365, 23, replace = TRUE)))>0}
sum(sapply(1:10000, birthdays))
}}}

{{{
과제: 어느 모임에서 생일이 같은 사람이 있을 확률이 99% 이상이 되려면 모임이 몇 명이어야 될까?
}}}

-----
<<Navigation(siblings,1)>>