Binomial self report

.gitignore

@@ -26,13 +26,16 @@ Manuscript/
 1_IGT_PP/Code/R/2_Plotting/group_PPCs_joint_indMu.Rmd
 1_IGT_PP/Code/Stan/simple_IRT.stan
 1_IGT_PP/Data/2_Fitted/IRT/
+1_IGT_PP/Figs_Tables/self_report/
 1_IGT_PP/Code/R/1_Analyses/Full_PANAS_simple_IRT.Rmd
 1_IGT_PP/Code/R/1_Analyses/NA_PANAS_simple_IRT.Rmd
 1_IGT_PP/Code/R/1_Analyses/PA_PANAS_simple_IRT.Rmd
 1_IGT_PP/Code/R/2_Plotting/posterior_predictive_checks_joint_indMu.Rmd
 1_IGT_PP/Code/R/1_Analyses/fitting_orl_joint_indMu.Rmd
 1_IGT_PP/Code/Stan/self_report_mus.stan
 1_IGT_PP/Code/R/1_Analyses/selfreport_indMu.RMD
+1_IGT_PP/Code/R/1_Analyses/selfreport_1mu.RMD
+1_IGT_PP/Code/R/1_Analyses/self_report_1mu.stan
 
 
 # 

1_IGT_PP/Code/R/0_Preprocess/Preprocess1.R → 1_IGT_PP/Code/R/0_Preprocess/IGT_Preprocess.R

@@ -104,7 +104,8 @@ AB <- mutate(AB, stim = recode(cardname, 'A'= 1, 'B' = 2, 'C' = 3, 'D' = 4), tim
 
 
 BD$cardname
-BD <- mutate(BD, stim = recode(cardname, 'A'= 3, 'B' = 1, 'C' = 4, 'D' = 2), time = SessionDate)
+# BD <- mutate(BD, stim = recode(cardname, 'A'= 3, 'B' = 1, 'C' = 4, 'D' = 2), time = SessionDate)
+BD <- mutate(BD, stim = recode(cardname, 'A'= 1, 'B' = 2, 'C' = 3, 'D' = 4), time = SessionDate)
 
 
 
@@ -267,10 +268,10 @@ unique(Sess2_IGT$ExperimentName)
 # Save IGT dataframes as .csv files
 ############################################################################################################
 
-Sess1_IGT_dat = here("Data", "1_Preprocessed", "Sess1_IGT.csv")
+Sess1_IGT_dat = here("1_IGT_PP", "Data", "1_Preprocessed", "Sess1_IGT.csv")
 write.csv(Sess1_IGT,Sess1_IGT_dat)
 
-Sess2_IGT_dat = here("Data", "1_Preprocessed", "Sess2_IGT.csv")
+Sess2_IGT_dat = here("1_IGT_PP", "Data", "1_Preprocessed", "Sess2_IGT.csv")
 write.csv(Sess2_IGT,Sess2_IGT_dat)
 
 
@@ -337,7 +338,7 @@ for (session in 1:2) {
   )
 
   # save Sess 1 or 2 lists as .rds
-  saveRDS(dataList, file = here("Data", "1_Preprocessed", paste0("Sess", session, ".rds")))
+  saveRDS(dataList, file = here("1_IGT_PP", "Data", "1_Preprocessed", paste0("Sess", session, ".rds")))
 }
 
 # I'm here, this is old code for preprocessing the raw data for running a joint ORL model
@@ -378,7 +379,8 @@ maxTrials <- max(Tsubj)
 # Behavioral data arrays
 # this creates 3 identical arrays of [1:50, 1:100, 1:2]; the array is filled with -1 values
 # array names are RLmatrix, SRLmatrix, choice
-RLmatrix <- SRLmatrix <- choice <- card <- array(-1, c(numSubjs, maxTrials,2))
+RLmatrix <- SRLmatrix <- choice <- card <- IDs <- array(-1, c(numSubjs, maxTrials,2))
+
 
 # Loop through and format into 3 arrays: choice, RLmatrix, SRLmatrix
 for (i in 1:numSubjs) {
@@ -392,11 +394,14 @@ for (i in 1:numSubjs) {
       card[i,,s] <- tmp_dat$stim
       RLmatrix[i,,s]  <- tmp_dat$rewlos
       SRLmatrix[i,,s] <- tmp_dat$Srewlos
+      IDs[i,,s] <- tmp_dat$Subject
     }
   }
 }
 
-
+# Check that how the for loop above loops through the subject IDs is in the same order as subjList
+(condensedIDs = apply(IDs, c(1,3), mean)) # Visually check they are the same (-1 means missing data on that session)
+ifelse(mean(condensedIDs[,1] == subjList) == 1, "Equal", "Not Equal") # Check if IDs in for loop equal IDs in subjList
 
 # Put in stan-ready list
 stan_dat <- list(
@@ -412,5 +417,5 @@ stan_dat <- list(
   subjIDs = subjList 
 )     
 
-saveRDS(stan_dat, file = here("1_IGT_PP", "Data", "1_Preprocessed", "stan_ready_ORL_joint_retest.rds"))
+saveRDS(stan_dat, file = here("1_IGT_PP", "Data", "1_Preprocessed", "stan_ready_ORL_IGT.rds"))
 

1_IGT_PP/Code/R/0_Preprocess/SelfReport_Preprocess.R

@@ -0,0 +1,176 @@
+
+
+
+library(dplyr)
+library(tidyverse) 
+library(haven)
+library(lubridate)
+library(readxl)
+library(here)
+
+
+
+# -----------------------------------------------
+# MEASURES
+# * Demographics
+#   + ID, sex, race, Ethnicity
+# 
+# * Behavioral Activation Scale
+#   + bastot    = ?
+#     + basdrive  = Drive subscale: measures persistent pursuit of goals
+#     + basfunsk  = Fun Seeking subscale: measures desire for and approach towards new rewards
+#     + basrewres = Reward Responsiveness subscale: measures positive responses to the anticipation or consummation of rewards
+#   
+# * Behavioral Inhibition Scale: measures sensitivity to negatively-valenced events
+#   + bis
+#   
+# * Positive Affect/Negative Affect Schedules: measures state-level mood
+#   + panas_pa = Positive Affect
+#   + panas_na = Negative Affect
+#   
+# * Mood and Anxiety Symptom Questionnaire (MASQ): assesses internalizing symptoms experienced in past week
+#   + masqGDA = General Distress Anxiety subscale: measures general anxious mood
+#   + masqAA  = Anxious Arousal subscale: measures somatic hyperarousal
+#   + masqGDD = General Distress Depression subscale: measures general depressed mood
+#   + masqAD  = Anhedonic Depression subscale: measures low positive affect
+#   
+# * Snaith-Hamilton Pleasure Scale: measures ability to experience pleasure in last few days
+#   + shaps_tot
+#   
+# * Patient-Reported Outcomes Measurement Information System: assesses symptoms of depression experienced in past week
+#   + prdep_tot
+
+
+
+# -----------------------------------------------
+## IMPORT DATA
+session1 = read_sav(here("1_IGT_PP", "Data", "0_Raw", "MergedQuest_3.21.16-Session1.sav"),
+                    col_select = c("ID", "sex", "race", "Ethnicity",
+                                   "bastot", "basdrive", "basfunsk",
+                                   "basrewres", "bis", "panas_pa",
+                                   "panas_na", "masqGDA", "masqAA", "masqGDD",
+                                   # "masqAD", # Removed
+                                   "shaps_tot", "prdep_tot")) %>% 
+  mutate(session = 1) %>%                     # Create session variable
+  filter(ID >= 2049,                          # Subset PP participants
+         ID != 2059,                          # Remove participant that played old IFT
+         ID != 2083)                          # Remove participant that had no task data
+
+session2 = read_sav(here("1_IGT_PP", "Data", "0_Raw", "MergedQuest_3.21.16-Session2.sav"),
+                    col_select = c("ID", "bastot", "basdrive", "basfunsk",
+                                   "basrewres", "bis", "panas_pa",
+                                   "panas_na", "shaps_tot")) %>% 
+  mutate(session = 2) %>%                     # Create session variable
+  filter(ID >= 2049,                          # Subset PP participants
+         ID != 2059,                          # Remove participant that played old IFT
+         ID != 2083,                          # Remove participant that had no task data
+         !ID %in% c(2057, 2063, 2064,         # Remove participants that didn't do session 2
+                    2067, 2086, 2090,
+                    2093, 2094, 2096, 2098))
+
+both_sessions = bind_rows(session1, session2) %>% 
+  select(ID, session,
+         # Session 1 only
+         prdep_tot, masqGDA, masqAA, masqGDD,
+         # masqAD, # Removed
+         # Both sessions
+         bastot, basdrive, basfunsk, basrewres,
+         bis, panas_pa, panas_na, shaps_tot,
+         # Remove demographics
+         -c(sex, race, Ethnicity)) %>% 
+  pivot_longer(cols = -c("ID", "session"), values_to = "score", names_to = "scale")
+
+
+
+# -----------------------------------------------
+# PREP DATA
+IDs = unique(both_sessions$ID)
+N = length(IDs)
+
+# Make dataframe describing each scale
+scales = data.frame(scale = c(# Session 1 only
+                              "prdep_tot", "masqGDA", "masqAA", "masqGDD",
+                              # "masqAD", # Removed
+                              # Both Sessions
+                              "bastot", "basdrive", "basfunsk", "basrewres",
+                              "bis", "panas_pa", "panas_na", "shaps_tot"),
+                    n_items = c(28, 11, 17, 12, 13, 4, 4, 5, 7, 10, 10, 14),
+                    min_item_score = c(0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0),
+                    max_item_score = c(4, 5, 5, 5, 4, 4, 4, 4, 4, 5, 5, 1),
+                    n_sessions = c(1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2)) %>% 
+  mutate(min = n_items * min_item_score,
+         max = n_items * max_item_score,
+         scaled_max = max - min)
+
+# Make scaled scores and check out which scores are below the minimums - for those scores, we'll truncate them to 0
+full_data = both_sessions %>% 
+  left_join(scales) %>% 
+  mutate(scaled_score = case_when(score-min < 0 ~ 0, is.numeric(score) ~ as.integer(round(score - min))))
+
+# To check for range problems, uncomment & run below
+# full_data %>%
+#   mutate(raw_pbx = case_when(score > max ~ "too high", score < min ~ "too low", T ~ ""),
+#          scaled_pbx = case_when(scaled_score > scaled_max ~ "too high", scaled_score < 0 ~ "too low", T ~ "")) %>%
+#   View(title = "Check_Pbxs")
+
+
+# Create lists to store data in for each scale
+participants = list()
+scores = list()
+missingness = list()
+maxes = list()
+
+for(l in 1:length(scales$scale)){
+  # Create matrices for the current scale in each list
+  participants[[scales$scale[l]]] = matrix(-999, nrow = N, ncol = scales$n_sessions[l])
+  scores[[scales$scale[l]]] = matrix(data = -999, nrow = N, ncol = scales$n_sessions[l])
+  missingness[[scales$scale[l]]] = matrix(data = 1, nrow = N, ncol = scales$n_sessions[l])
+
+  maxes[[scales$scale[l]]] = filter(scales, scale == scales$scale[l])[,"scaled_max"] # Save max score for current scale
+
+  cur_scale = filter(full_data, scale == scales$scale[l]) # Subset scores from scale of interest
+
+  for(s in 1:scales$n_sessions[l]){
+    cur_session = filter(cur_scale, session == s) # Subset scores from session of interest
+
+    for(i in 1:N){
+      cur_obs = filter(cur_session, ID == IDs[i]) # Subset score from person of interest
+      participants[[scales$scale[l]]][i,s] = IDs[i]
+
+      if(length(cur_obs$scaled_score)){ # Check if data is missing, if not, then get data
+        scores[[scales$scale[l]]][i,s] = cur_obs$scaled_score # Save score in element of matrix within the scale's score matrix
+        missingness[[scales$scale[l]]][i,s] = 0 # Save data as not missing in element of matrix within the scale's missingness matrix
+      }
+    }
+  }
+}
+
+
+stan_datas = list()
+
+for(scale in unique(scales$scale)){
+  stan_datas[[scale]] = list(IDs = participants[[scale]],
+                             N = nrow(scores[[scale]]),       # Number of participants for current scale
+                             S = ncol(scores[[scale]]),       # Number of sessions for current scale
+                             missing = missingness[[scale]],  # Matrix of missingness for current scale
+                             score = scores[[scale]],         # Matrix of scores for current scale
+                             M = maxes[[scale]])              # Max score for current scale
+}
+
+
+saveRDS(stan_datas, file = here("1_IGT_PP", "Data", "1_Preprocessed", "stan_ready_binomial_selfreport.rds"))
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+