# Roger Palay copyright 2016-02-20
# Saline, MI 48176
#
ci_2unknown <- function (
  s_one, n_one, x_one,
  s_two, n_two, x_two,
  cl=0.95)
{
    # try to avoid some common errors
    if( (cl <=0) | (cl>=1) )
      {return("Confidence interval must be strictly between 0.0 and 1")
      }
    if( (s_one <= 0)  | (s_two <= 0) )
      {return("Sample standard deviation must be positive")}
    if( (as.integer(n_one) != n_one ) |
        (as.integer(n_two) != n_two ) )
      {return("Sample size must be a whole number")}

    alpha <- 1 - cl
    alpha_div_2 <- alpha / 2

    d_one <- s_one^2/n_one
    d_two <- s_two^2/n_two
    std_err <- sqrt( d_one+d_two )

    if( n_one < n_two)
       { df_simple <- n_one - 1} else
       { df_simple <- n_two - 1}

    t_simple <- qt(alpha_div_2,
                   df_simple,lower.tail=FALSE)

    #compute the degrees of freedom
    d_one <- s_one^2/n_one
    d_two <- s_two^2/n_two

    df <-(d_one+d_two)^2/
         ( d_one^2/(n_one-1)+d_two^2/(n_two-1))

    t_full <- qt(alpha_div_2,
                 df,lower.tail=FALSE)

    diff <- x_one - x_two
    # our point estimate

    simple_moe <- t_simple*std_err

    simple_low <- diff - simple_moe
    simple_high<- diff + simple_moe

    full_moe <- t_full*std_err

    full_low <- diff - full_moe
    full_high <- diff + full_moe

    result <- c(full_low, full_high, full_moe,
                simple_low, simple_high, simple_moe,
                std_err,  alpha_div_2,  df, t_full,
                df_simple, t_simple, diff )

    names(result)<-c("Full Low", "Full High", "Full MOE",
                     "Simp Low", "Simp High", "Simp MOE",
                     "St. Error", "alpha/2","DF calc",
                     "t calc", "DF Simp", "t Simp", "Diff")
    return( result )
}