Make summarize estimator group correctly for robustness checks.

[ml_measurement_error_public.git] / civil_comments / design_example.R

diff --git a/civil_comments/design_example.R b/civil_comments/design_example.R
index 5991334c3cbbdc69649da5bdce4d835a8c95eced..49076880fb5f68ff362dc6730fdec8ac99f542b1 100644 (file)
--- a/civil_comments/design_example.R
+++ b/civil_comments/design_example.R
@@ -1,126 +1,10 @@
-library(data.table)
-library(MASS)
-
-scores <- fread("perspective_scores.csv")
-scores <- scores[,id:=as.character(id)]
-
-df <- fread("all_data.csv")
-
-# only use the data that has identity annotations
-df <- df[identity_annotator_count > 0]
-
-(df[!(df$id %in% scores$id)])
-
-df <- df[scores,on='id',nomatch=NULL]
-
+set.seed(1111)
+source('load_perspective_data.R')

 ## how accurate are the classifiers?
 
 ## the API claims that these scores are "probabilities"
 ## say we care about the model of the classification, not the probability
 
 ## how accurate are the classifiers?
 
 ## the API claims that these scores are "probabilities"
 ## say we care about the model of the classification, not the probability
 

-F1 <- function(y, predictions){
-    tp <- sum( (predictions == y) & (predictions==1))
-    fn <- sum( (predictions != y) & (predictions!=1))
-    fp <- sum( (predictions != y) & (predictions==1))
-    precision <- tp / (tp + fp)
-    recall <- tp / (tp + fn)
-    return (2 * precision * recall ) / (precision + recall)
-}
-
-df[, ":="(identity_attack_pred = identity_attack_prob >=0.5,
-          insult_pred = insult_prob >= 0.5,
-          profanity_pred = profanity_prob >= 0.5,
-          severe_toxicity_pred = severe_toxicity_prob >= 0.5,
-          threat_pred = threat_prob >= 0.5,
-          toxicity_pred = toxicity_prob >= 0.5,
-          identity_attack_coded = identity_attack >= 0.5,
-          insult_coded = insult >= 0.5,
-          profanity_coded = obscene >= 0.5,
-          severe_toxicity_coded = severe_toxicity >= 0.5,
-          threat_coded = threat >= 0.5,
-          toxicity_coded = toxicity >= 0.5
-          )]
-
-
-
-## toxicity is about 93% accurate, with an f1 of 0.8
-## identity_attack has high accuracy 97%, but an unfortunant f1 of 0.5.
-## threat has high accuracy 99%, but a really bad looking f1 of 0.48.
-accuracies <- df[,.(identity_attack_acc = mean(identity_attack_pred == identity_attack_coded),
-                    insult_pred_acc = mean(insult_pred == insult_coded),
-                    profanity_acc = mean(profanity_pred == profanity_coded),
-                    severe_toxicity_acc = mean(severe_toxicity_pred == severe_toxicity_coded),
-                    theat_acc = mean(threat_pred == threat_coded),
-                    toxicity_acc = mean(toxicity_pred == toxicity_coded))]
-
-f1s <- df[,.(identity_attack_f1 = F1(identity_attack_coded,identity_attack_pred),
-                    insult_f1 = F1(insult_coded,insult_pred),
-                    profanity_f1 = F1(profanity_coded,profanity_pred),
-                    severe_toxicity_f1 = F1(severe_toxicity_coded,severe_toxicity_pred),
-                    theat_f1 = F1(threat_coded,threat_pred),
-                    toxicity_f1 = F1(toxicity_coded,toxicity_pred))]
-
-positive_cases <- df[,.(identity_attacks = sum(identity_attack_coded),
-                        insults = sum(insult_coded),
-                        profanities = sum(profanity_coded),
-                        severe_toxic_comments = sum(severe_toxicity_coded),
-                        threats = sum(threat_coded),
-                        toxic_comments = sum(toxicity_coded))]
-
-## there are 50,000 toxic comments, 13000 identity attacks, 30000 insults, 3000 profanities, 8 severe toxic, and 1000 threats.
-
-proportions_cases <- df[,.(prop_identity = mean(identity_attack_coded),
-                           prop_insults = mean(insult_coded),
-                           prop_profanity = mean(profanity_coded),
-                           prop_severe = mean(severe_toxicity_coded),
-                           prop_threats = mean(threat_coded),
-                           prop_toxic = mean(toxicity_coded))]
-
-## at 11% of comments, "toxicity" seems not so badly skewed. Try toxicity first, and if it doesn't work out try insults.
-
-## now look for an example where differential error affects an identity, or a reaction.
-df <- df[,":="(identity_error = identity_attack_coded - identity_attack_pred,
-               insult_error = insult_coded - insult_pred,
-               profanity_error = profanity_coded - profanity_pred,
-               severe_toxic_error = severe_toxicity_coded - severe_toxicity_pred,
-               threat_error = threat_coded - threat_pred,
-               toxicity_error = toxicity_coded - toxicity_pred)]
-
-## what's correlated with toxicity_error ?
-df <- df[,approved := rating == "approved"]
-
-cortab <- cor(df[,.(toxicity_error,
-                    identity_error,
-                    toxicity_coded,
-                    funny,
-                    approved,
-                    sad,
-                    wow,
-                    likes,
-                    disagree,
-                    male,
-                    female,
-                    transgender,
-                    other_gender,
-                    heterosexual,
-                    bisexual,
-                    other_sexual_orientation,
-                    christian,
-                    jewish,
-                    hindu,
-                    buddhist,
-                    atheist,
-                    other_religion,
-                    black,
-                    white,
-                    asian,
-                    latino,
-                    other_race_or_ethnicity,
-                    physical_disability,
-                    intellectual_or_learning_disability,
-                    psychiatric_or_mental_illness,
-                    other_disability)])
-

 ## toxicity error is weakly correlated pearson's R = 0.1 with both "white" and "black".
 
 ## compare regressions with "white" or "black" as the outcome and "toxicity_coded" or "toxicity_pred" as a predictor.
 ## toxicity error is weakly correlated pearson's R = 0.1 with both "white" and "black".
 
 ## compare regressions with "white" or "black" as the outcome and "toxicity_coded" or "toxicity_pred" as a predictor.


@@ -134,14 +18,46 @@ cortab['toxicity_coded',]
 cortab['identity_error',]
 cortab['white',]
 
 cortab['identity_error',]
 cortab['white',]
 

-glm(white ~ toxicity_coded + psychiatric_or_mental_illness, data = df, family=binomial(link='logit'))
+## here's a simple example, is P(white | toxic and mentally ill) > P(white | toxic or mentally ill). Are people who discuss their mental illness in a toxic way more likely to be white compared to those who just talk about their mental illness or are toxic? 
+summary(glm(white ~ toxicity_coded*psychiatric_or_mental_illness, data = df, family=binomial(link='logit')))
+
+summary(glm(white ~ toxicity_pred*psychiatric_or_mental_illness, data = df, family=binomial(link='logit')))

 
 

-glm(white ~ toxicity_pred + psychiatric_or_mental_illness, data = df, family=binomial(link='logit'))
+summary(glm(white ~ toxicity_coded*male, data = df, family=binomial(link='logit')))

 
 

-m1 <- glm.nb(funny ~ (male + female + transgender + other_gender + heterosexual + bisexual + other_sexual_orientation + christian + jewish + hindu + buddhist + atheist + other_religion + asian + latino + other_race_or_ethnicity + physical_disability + intellectual_or_learning_disability + white + black + psychiatric_or_mental_illness)*toxicity_coded, data = df)
+summary(glm(white ~ toxicity_pred*male, data = df, family=binomial(link='logit')))
+
+summary(glm(toxicity_coded ~ white*psychiatric_or_mental_illness, data = df, family=binomial(link='logit')))
+
+summary(glm(toxicity_pred ~ white*psychiatric_or_mental_illness, data = df, family=binomial(link='logit')))
+
+
+## another simple enough example: is P(toxic | funny and white) > P(toxic | funny nand white)? Or, are funny comments more toxic when people disclose that they are white?
+
+summary(glm(toxicity_pred ~ funny*white, data=df, family=binomial(link='logit')))
+summary(glm(toxicity_coded ~ funny*white, data=df, family=binomial(link='logit')))
+
+source("../simulations/measerr_methods.R")
+                                                                                       
+saved_model_file <- "measerr_model_tox.eq.funny.cross.white.RDS"
+overwrite_model <- TRUE 
+
+# it works so far with a 20% and 15% sample. Smaller is better. let's try a 10% sample again. It didn't work out. We'll go forward with a 15% sample.
+df_measerr_method <- copy(df)[sample(1:.N, 0.05 * .N), toxicity_coded_1 := toxicity_coded]
+df_measerr_method <- df_measerr_method[,toxicity_coded := toxicity_coded_1]
+summary(glm(toxicity_coded ~ funny*white, data=df_measerr_method[!is.na(toxicity_coded)], family=binomial(link='logit')))
+
+if(!file.exists(saved_model_file) || (overwrite_model == TRUE)){
+    measerr_model <- measerr_mle_dv(df_measerr_method,toxicity_coded ~ funny*white,outcome_family=binomial(link='logit'), proxy_formula=toxicity_pred ~ toxicity_coded*funny*white)
+saveRDS(measerr_model, saved_model_file)
+} else {
+    measerr_model <- readRDS(saved_model_file)
+}

 
 

-m2 <- glm.nb(funny ~ (male + female + transgender + other_gender + heterosexual + bisexual + other_sexual_orientation + christian + jewish + hindu + buddhist + atheist + other_religion + asian + latino + other_race_or_ethnicity + physical_disability + intellectual_or_learning_disability + white + black + psychiatric_or_mental_illness)*toxicity_pred, data = df)
+inv_hessian <- solve(measerr_model$hessian)
+se <- diag(inv_hessian)

 
 

+lm2 <- glm.nb(funny ~ (male + female + transgender + other_gender + heterosexual + bisexual + other_sexual_orientation + christian + jewish + hindu + buddhist + atheist + other_religion + asian + latino + other_race_or_ethnicity + physical_disability + intellectual_or_learning_disability + white + black + psychiatric_or_mental_illness)*toxicity_pred, data = df)

 m3 <- glm.nb(funny ~ (male + female + transgender + other_gender + heterosexual + bisexual + other_sexual_orientation + christian + jewish + hindu + buddhist + atheist + other_religion + asian + latino + other_race_or_ethnicity + physical_disability + intellectual_or_learning_disability + white + black + psychiatric_or_mental_illness)*toxicity, data = df)
 
 
 m3 <- glm.nb(funny ~ (male + female + transgender + other_gender + heterosexual + bisexual + other_sexual_orientation + christian + jewish + hindu + buddhist + atheist + other_religion + asian + latino + other_race_or_ethnicity + physical_disability + intellectual_or_learning_disability + white + black + psychiatric_or_mental_illness)*toxicity, data = df)