3737#
3838# 1. Changed from one-sided test to two-sided test in all t-tests
3939
40+ # ###################################
41+ # Update by Tianchen Qian, 2019/3/21
42+ #
43+ # 1. Add alpha and alpha_se to the output.
44+
45+
4046
4147library(rootSolve ) # for solver function multiroot()
4248
@@ -426,7 +432,9 @@ SARA_primary_hypothesis_1 <- function(
426432 beta_se = as.numeric(result $ beta_se_ssa ),
427433 test_stat = as.numeric(result $ test_result_t $ test_stat ),
428434 critical_value = result $ test_result_t $ critical_value ,
429- p_value = as.numeric(result $ test_result_t $ p_value ))
435+ p_value = as.numeric(result $ test_result_t $ p_value ),
436+ alpha = as.numeric(result $ alpha_hat ),
437+ alpha_se = as.numeric(result $ alpha_se_ssa ))
430438 return (output )
431439}
432440
@@ -522,7 +530,9 @@ SARA_primary_hypothesis_2 <- function(
522530 beta_se = as.numeric(result $ beta_se_ssa ),
523531 test_stat = as.numeric(result $ test_result_t $ test_stat ),
524532 critical_value = result $ test_result_t $ critical_value ,
525- p_value = as.numeric(result $ test_result_t $ p_value ))
533+ p_value = as.numeric(result $ test_result_t $ p_value ),
534+ alpha = as.numeric(result $ alpha_hat ),
535+ alpha_se = as.numeric(result $ alpha_se_ssa ))
526536 return (output )
527537}
528538
@@ -591,14 +601,16 @@ SARA_exploratory_analysis <- function(
591601 avail_var = avail_var ,
592602 prob_treatment = prob_treatment ,
593603 significance_level = significance_level )
594- output <- list (beta = result $ beta_hat ,
595- beta_se = result $ beta_se_ssa ,
604+ output <- list (beta = as.numeric( result $ beta_hat ) ,
605+ beta_se = as.numeric( result $ beta_se_ssa ) ,
596606 test_stat_t = as.numeric(result $ test_result_t $ test_stat ),
597607 critical_value_t = result $ test_result_t $ critical_value ,
598608 p_value_t = as.numeric(result $ test_result_t $ p_value ),
599609 test_stat_f = result $ test_result_f $ test_stat ,
600610 critical_value_f = result $ test_result_f $ critical_value ,
601- p_value_f = result $ test_result_f $ p_value )
611+ p_value_f = result $ test_result_f $ p_value ,
612+ alpha = as.numeric(result $ alpha_hat ),
613+ alpha_se = as.numeric(result $ alpha_se_ssa ))
602614 return (output )
603615}
604616
@@ -704,14 +716,16 @@ SARA_exploratory_analysis_general_F_test <- function(
704716 p_value_f <- pf(test_stat_f , df1 = p1 , df2 = n - q - p , lower.tail = FALSE )
705717
706718
707- output <- list (beta = result $ beta_hat ,
708- beta_se = result $ beta_se_ssa ,
719+ output <- list (beta = as.numeric( result $ beta_hat ) ,
720+ beta_se = as.numeric( result $ beta_se_ssa ) ,
709721 test_stat_t = as.numeric(result $ test_result_t $ test_stat ),
710722 critical_value_t = result $ test_result_t $ critical_value ,
711723 p_value_t = as.numeric(result $ test_result_t $ p_value ),
712724 test_stat_f = test_stat_f ,
713725 critical_value_f = critical_value_f ,
714- p_value_f = p_value_f )
726+ p_value_f = p_value_f ,
727+ alpha = as.numeric(result $ alpha_hat ),
728+ alpha_se = as.numeric(result $ alpha_se_ssa ))
715729 return (output )
716730}
717731
@@ -791,6 +805,12 @@ if (0) {
791805 $ p_value
792806 [1 ] 1.022584e-14
793807
808+ $ alpha
809+ [1 ] - 0.4286552699 0.0001338947 0.0022651338
810+
811+ $ alpha_se
812+ [1 ] 0.0271394442 0.0292520053 0.0003041147
813+
794814 >
795815 > # primary hypothesis 2
796816 > SARA_primary_hypothesis_2(dta , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), survey_completion_var = " Y" )
@@ -809,6 +829,12 @@ if (0) {
809829 $ p_value
810830 [1 ] 0.7238992
811831
832+ $ alpha
833+ [1 ] - 0.1914996030 - 0.0078097155 0.0001918152
834+
835+ $ alpha_se
836+ [1 ] 0.0287537133 0.0294240517 0.0003142004
837+
812838 >
813839 > # exploratory analysis
814840 > SARA_exploratory_analysis(dta , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), moderator = " Y_lag1" )
@@ -838,6 +864,12 @@ if (0) {
838864 $ p_value_f
839865 [1 ] 1.81756e-23
840866
867+ $ alpha
868+ [1 ] - 0.522625830 0.121815023 0.002169604
869+
870+ $ alpha_se
871+ [1 ] 0.0449813807 0.0515054456 0.0003145783
872+
841873 > # this will give the same F-test result as SARA_exploratory_analysis()
842874 > SARA_exploratory_analysis_general_F_test(dta , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), moderator = " Y_lag1" , F_test_L = diag(2 ))
843875 $ beta
@@ -866,6 +898,12 @@ if (0) {
866898 $ p_value_f
867899 [1 ] 1.81756e-23
868900
901+ $ alpha
902+ [1 ] - 0.522625830 0.121815023 0.002169604
903+
904+ $ alpha_se
905+ [1 ] 0.0449813807 0.0515054456 0.0003145783
906+
869907 > SARA_exploratory_analysis_general_F_test(dta , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), moderator = " Y_lag1" , F_test_L = rep(1 , 2 ))
870908 $ beta
871909 Intercept Y_lag1
@@ -893,6 +931,12 @@ if (0) {
893931 $ p_value_f
894932 [1 ] 6.113323e-09
895933
934+ $ alpha
935+ [1 ] - 0.522625830 0.121815023 0.002169604
936+
937+ $ alpha_se
938+ [1 ] 0.0449813807 0.0515054456 0.0003145783
939+
896940 >
897941 >
898942 > # ## create fake availability indicator, and try the three analysis functions with availability ###
@@ -916,6 +960,12 @@ if (0) {
916960 $ p_value
917961 [1 ] 0.03407131
918962
963+ $ alpha
964+ [1 ] - 0.32417513 - 0.03184334 0.00195437
965+
966+ $ alpha_se
967+ [1 ] 0.0563985798 0.0640324410 0.0006506912
968+
919969 > SARA_primary_hypothesis_2(dta2 , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), survey_completion_var = " Y" , avail_var = " avail" )
920970 $ beta
921971 [1 ] 0.02120382
@@ -932,6 +982,12 @@ if (0) {
932982 $ p_value
933983 [1 ] 0.6068408
934984
985+ $ alpha
986+ [1 ] - 0.2131429183 0.0033599799 0.0000873171
987+
988+ $ alpha_se
989+ [1 ] 0.0485193353 0.0605234082 0.0008457411
990+
935991 > SARA_exploratory_analysis(dta2 , control_var = c(" Y_lag1" , " at_tapcount_lag1" ), moderator = " Y_lag1" , avail_var = " avail" )
936992 $ beta
937993 Intercept Y_lag1
@@ -958,4 +1014,10 @@ if (0) {
9581014
9591015 $ p_value_f
9601016 [1 ] 0.01252342
961- }
1017+
1018+ $ alpha
1019+ [1 ] - 0.352447929 0.006808953 0.001907070
1020+
1021+ $ alpha_se
1022+ [1 ] 0.0873681644 0.1024315444 0.0006620606
1023+ }
0 commit comments