임시

#1. 일표본 t검정
#data는 벡터형식
#pmean은 모집단평균
#양측검정
t.test(data,mu=pmean)
#단측 우측꼬리검정
t.test(data,mu=pmean,alternative="greater")
#단측 좌측꼬리검정
t.test(data,mu=pmean,alternative="less")

#2. 독립표본 t검정(등분산 가정)
#data1과 data2는 벡터형식
#양측검정
t.test(data1,data2,var.equal=TRUE)
#단측 우측꼬리검정
t.test(data1,data2,var.equal=TRUE,alternative="greater")
#단측 좌측꼬리검정
t.test(data1,data2,var.equal=TRUE,alternative="less")

#3. 독립표본 t검정(이분산 가정)
#data1과 data2는 벡터형식
#양측검정
t.test(data1,data2)
#단측 우측꼬리검정
t.test(data1,data2,alternative="greater")
#단측 좌측꼬리검정
t.test(data1,data2,alternative="less")

#4. 대응표본 t검정
#data1과 data2는 벡터형식
#양측검정
t.test(data1,data2,paired=TRUE)
#단측 우측꼬리검정
t.test(data1,data2,paired=TRUE,alternative="greater")
#단측 좌측꼬리검정
t.test(data1,data2,paired=TRUE,alternative="less")

#5. 윌콕슨 순위합 검정(독립표본, 비모수)
#data1과 data2는 벡터형식
#양측검정
wilcox.test(data1,data2)
#단측 우측꼬리검정
wilcox.test(data1,data2,alternative="greater")
#단측 좌측꼬리검정
wilcox.test(data1,data2,alternative="less")

#6. 윌콕슨 부호순위 검정(대응표본, 비모수)
#data1과 data2는 벡터형식
#양측검정
wilcox.test(data1,data2,paired=TRUE)
#단측 우측꼬리검정
wilcox.test(data1,data2,paired=TRUE,alternative="greater")
#단측 좌측꼬리검정
wilcox.test(data1,data2,paired=TRUE,alternative="less")

#R로 크루스칼-왈리스 검정해보기(kruskal.test() 함수)
> A <- c(257, 205, 206, 231, 190, 214, 228, 203)
> B <- c(201, 164, 197, 185)
> C <- c(248, 165, 187, 220, 212, 215, 281)
> D <- c(202, 276, 207, 204, 230, 227)
> kruskal.test(list(A,B,C,D)) #결과는 
Kruskal-Wallis rank sum test
data:  list(A, B, C, D)
Kruskal-Wallis chi-squared = 6.839, df = 3, p-value = 0.07721

# 123, 105, 117, 117,109, 118, 122 의 중위수는 118 이라고 주장하는 경우 부호검정
theta <- 118 #주의: 표본의 중위수가 아닌 estimated median value
data <- c(123, 105, 117, 117,109, 118, 122)

# 부호검정
library(BSDA)
SIGN.test(data, md=theta, alternative = "two.sided")

# 이항분포를 이용한 검정
S <- sum(data > theta)    # theta 보다 큰 갯수
n <- sum(data != theta)    # n : theta 가 아닌 갯수. 즉 유효한 샘플수
s <- min(n-S,S)                 # s : + - 중 작은 부호의 갯수
pbinom(s, n, 0.5, lower.tail = TRUE) * 2   # n 개중에 s 개를 뽀는 이항분포 확률

library(MASS) 
library(caret) 

# createDataPartition 
set.seed(123) 
train <- createDataPartition(y=Boston$medv, p=0.7, list=F) 
head(train) 

Boston.train <- Boston[train,] 
Boston.test <- Boston[-train,] 

#glmnet elasticnet : alpha = 0.1 lambda = 0.08747614.
set.seed(123)
Boston.gnet <- train(form=medv ~ ., data=Boston.train, method="glmnet",  
                   trControl=trainControl(method="cv",number=10), tuneLength=10)  

#rf 오래걸림 : mtry = 5
Boston.rf <- train(form=medv ~ ., data=Boston.train, method="rf",  
                   trControl=trainControl(method="cv",number=10), tuneLength=10)

#svm : sigma = 0.04434552 and C = 32.
Boston.svm <- train(form=medv ~ ., data=Boston.train, method="svmRadialSigma",  
                   trControl=trainControl(method="cv",number=10), tuneLength=10)

# predict
gnet.pred <- predict(Boston.gnet,Boston.test)     
postResample(pred=gnet.pred, obs=Boston.test$medv) 

rf.pred <- predict(Boston.rf,Boston.test)     
postResample(pred=rf.pred, obs=Boston.test$medv) 

svm.pred <- predict(Boston.svm,Boston.test)     
postResample(pred=svm.pred, obs=Boston.test$medv) 


#comparison
models <- list(gnet=Boston.gnet, rf=Boston.rf, svm=Boston.svm) 
summary(resamples(models)) 
summary(resamples(models), metric="RMSE") 
summary(diff(resamples(models), metric="RMSE"))

#2
m <- matrix(c(15,60,25,25,69,6,10,77,13), ncol=3) 
dimnames(m) <- list(score=c("1.5-2.5","2.5-3.5","3.5-4.5"),
                    subject=c("사회과학","자연과학","공학")) 
m
addmargins(m) #합계 확인
Xsq <- chisq.test(m)
Xsq$expected
#3
chisq.test(m)

#연속형 확률분포 (Continuous probability distribution)에는
#정규분포: norm(), 균등분포: unif(), 지수분포: exp(), t-분포: t(), F-분포: f(), 카이제곱분포: chisq()

#누적 t분포 확률 값 계산 : pt(q, df, lower.tail=TRUE/FALSE) 
pt(q=1, df=1, lower.tail = TRUE) #0.75 pt는 누적분포함수, dt는 밀도함수
# t분포 분위수 함수 값 계산 : qt(p, df, lower.tail=TRUE/FALSE)
qt(p=0.75, df=1, lower.tail = TRUE) #1
#자유도가 1인 t분포에서,  50개의 실수값을 임의추출
rt(50, df=1)

#같은 흐름으로 for 정규분포
pnorm(q, mean=0, sd=1, lower.tail=TRUE/FALSE) #cf. dnorm(x, mean=0, sd=1)
qnorm(p, mean=0, sd=1, lower.tail=TRUE/FALSE)
rnorm(n, mean=0, sd=1)

#for F분포
pf(df1, df2, lower.tail=TRUE/FALSE) #cf. df(df1, df2)
qf(df1, df2, lower.tail=TRUE/FALSE)
rf(n, df1, df2)

library(TTR)
library(forecast)
library(tseries)
king <- scan("https://robjhyndman.com/tsdldata/misc/kings.dat", skip=3)
king.ts <- ts(king)

train <- subset(king.ts, end=length(king.ts)-4)
test  <- subset(king.ts, start=length(king.ts)-3)

#총42년의 데이터중 약10%인 마지막 4년을 테스트데이터로 사용
king.ts
train <- subset(king.ts, end=length(king.ts)-4)
test  <- subset(king.ts, start=length(king.ts)-3)

# 1. 정상성 체크
#정상성 평가 - 비정상
adf.test(train)
# 필요한 차분수 확인 - 1차 차분
ndiffs(train)
#1차 차분 데이터의 정상성 평가 - 정상
adf.test(diff(train))

#1차 차분 pq 수치 확인 -> pq 가 3,1
par(mfrow=c(1,2))
acf(diff(train));pacf(diff(train))

# 2. 시계열 예측 모델 제시
# pdq  가 3,1,1 이므로 다음 3개중 하나 선택
m1 <- arima(king.ts, order=c(3,1,0))
m2 <- arima(king.ts, order=c(0,1,1))
m3 <- arima(king.ts, order=c(3,1,1))

 

# 3. 한가지 모델을 최종 선택
# RMSE 관점에서 3,1,1 모델이 가장 우수
# 단, 간명도 원칙에 따라 차이가 크지 않을 경우 가장 간단한 모수를 설정
library(dplyr)
rbind(cbind(m='m1',accuracy(m1)),
      cbind(m='m2',accuracy(m2)),
      cbind(m='m3',accuracy(m3))) %>%
  data.frame() %>% 
  select(m,RMSE) %>% 
  arrange(RMSE)

#정규성 평가
qqnorm(m3$residuals, pch=21, col="black", 
       bg="gold", main="Q-Q Plot of Residuals") 

qqline(m3$residuals, col="royalblue", lwd=2) 
# H0 : 데이터들이 독립적으로 분포함. 
# p-value = 0.9349 그러므로 독립
Box.test(m3$residuals, type="Ljung-Box") 

 

# 4. 최종 예측 후 실제 결과와 비교 평가
# 예측
m3.pred <- forecast(m3,4)
plot(m3.pred)
plot(m3.pred, 
     col="darkgreen", 
     lwd=2, 
     flty=1,flwd=2, 
     fcol="royalblue", shadecols=c("mistyrose","salmon"), 
     xlab="y", 
     ylab="x", 
     main="Forecast" 
)

#x구성
library(dplyr)
set.seed(3312)
round(rnorm(40, mean = 10, sd = 4)) %>%
  data.frame(errors=.) %>%
  mutate(lot=row_number()) %>%
  mutate(errors=ifelse(errors<0,0,errors)) %>% 
  select(2,1) -> x
#LCL, CL, UCL 구성
x
summary(x)
plot(x)
sum(x$errors) 
pbar <- sum(x$errors) / (200 * 40) #에러율의 평균
sigma <- sqrt(pbar *(1-pbar) / 200) #확률의 표준편차
CL <- pbar
UCL <- pbar+3*sigma
LCL <- pbar-3*sigma
#시각화
library(ggplot2)
x %>% mutate(r=(errors/200))  %>%
  ggplot(aes(lot, r)) +
  geom_line(size=1, color = 'blue') +
  geom_line(aes(y = UCL),col = 'red') + # 관리상한선
  geom_line(aes(y = CL),col = 'red') + # 관리중심선
  geom_line(aes(y = LCL),col = 'red') # 관리하한선

#R
# Install the required package
install.packages("fastDummies")
  
# Load the library
library(fastDummies)
  
# Using PlantGrowth dataset
data <- PlantGrowth
  
# Create dummy variable
data <- dummy_cols(data, 
                   select_columns = "group")

#Python
X = pd.get_dummies(Data = X, columns=['group']

n <- 100 # 표본공간
x <- 80 # 성공횟수
p <- x/n   # 성공확률
prop.test(x,n,p,alternative="two.sided")