# encording with Shift-JIS

#
# Chapter 6. Test
#

## C:\EconDat\6_test
setwd("C:/EconDat/6_test")
getwd()


#
# 例6.3 大学生の平均アルバイト時給の仮説検定のための検定統計量
#

# 学生生活アンケート2018
df_survey <- read.csv("../data/student_survey2018/student_survey2018_data.csv", header=TRUE)

#データの構造を知る
str(df_survey)

id <- df_survey$id

#Q45	アルバイトの時給はいくらですか？(単位：円)
wage <- df_survey$Q45

# Q2	あなたの性別を教えて下さい。
sex <- df_survey$Q2

#データフレームの作成
df_wage <- data.frame(id, wage, sex)
length(df_wage$wage)

# 欠損値の削除
df_wage <- na.omit(df_wage)
length(df_wage$wage)


# 平均賃金
wage_mean <- mean(df_wage$wage)
wage_sd <- sd(df_wage$wage)
wage_n <- length(df_wage$wage)

#標本平均の標準誤差
wage_mean_se <- wage_sd/sqrt(wage_n)
wage_mean_se

# t値
wage_mean_t <- (wage_mean-950)/wage_mean_se
wage_mean_t


#
# 例6.7 t.test関数を用いた大学生の平均アルバイト時給の両側検定
#

## 両側検定
## 有意水準5%
t.test(wage, mu=950, alternative="two.sided")

## 有意水準1%
t.test(wage, mu=950, alternative="two.sided", conf.level=0.99)

#
# 例6.8 大学生の平均アルバイト時給の両側検定のt分布の臨界値
#

# 自由度の計算
df = wage_n -1
df

# 下側の臨界値
t_lower <- qt(0.025,df)
t_lower

# 上側の臨界値
t_upper <- qt(0.975,df)
t_upper

# 95%信頼区間

wage_mean_lower <- wage_mean+t_lower*wage_mean_se
wage_mean_lower

wage_mean_upper <- wage_mean-t_upper*wage_mean_se
wage_mean_upper

#
# 付論 Appendix 
# A6 6章 仮説検定の考え方・付論
#

#
# 例6.9：大学生の平均アルバイト時給の仮説
#

# 臨界値
t_upper <- qt(0.95,df)
t_upper

#片側検定
t.test(wage, mu=950, alternative="greater")

#
# 例6.10: 大学生アルバイト時給の男女差
#

# 男女別賃金
## 平均
wage_mean_sex <-tapply(df_wage$wage, df_wage$sex, mean)
wage_mean_sex

## 標準偏差
wage_sd_sex <-tapply(df_wage$wage, df_wage$sex, sd)
wage_sd_sex

## 標本サイズ
wage_n_sex <-tapply(df_wage$wage, df_wage$sex, length)
wage_n_sex

#男女別に記述統計量を保存
#平均
wage_male_mean <- wage_mean_sex["男性"]
wage_male_mean

wage_female_mean <- wage_mean_sex["女性"]
wage_female_mean

# 標準偏差
wage_male_sd <- wage_sd_sex["男性"]
wage_male_sd

wage_female_sd <- wage_sd_sex["女性"]
wage_female_sd

#標本サイズ
wage_male_n <- wage_n_sex["男性"]
wage_male_n

wage_female_n <- wage_n_sex["女性"]
wage_female_n

#標本平均の標準誤差
wage_male_mean_se <- wage_male_sd/sqrt(wage_male_n)
wage_male_mean_se

wage_female_mean_se <- wage_female_sd/sqrt(wage_female_n)
wage_female_mean_se


#
# 例6.11: 大学生アルバイト時給の男女差の仮説検定
#

#　平均値の差の検定
## 平均値の差の標準誤差
t_diff_se <- sqrt(wage_male_mean_se^2 +wage_female_mean_se^2)
## 検定統計量
t_diff <- (wage_male_mean - wage_female_mean)/t_diff_se
names(t_diff) <- "平均差"
t_diff

## 臨界値
qnorm(0.025)
qnorm(0.975)

#
# 例6.12：t.test関数を用いた大学生アルバイト時給の男女差の仮説検定
#

#　t.test関数を用いた平均値の差の検定
## 男性のアルバイト時給
wage_male <- df_wage$wage[df_wage$sex=="男性"]
## 女性のアルバイト時給
wage_female <- df_wage$wage[df_wage$sex=="女性"]

## 平均値の差の検定
t.test(wage_male, wage_female, mu=0,  var.equal=FALSE, alternative="two.sided") 

## 平均値の差の検定(Welchの検定)での自由度の計算
numerator = (wage_male_sd^2/wage_male_n +wage_female_sd^2/wage_female_n )^2
denominator =wage_male_sd^4/(wage_male_n^2*(wage_male_n-1)) +wage_female_sd^4/(wage_female_n^2*(wage_female_n-1))

t_diff_df = numerator/denominator
t_diff_df
