R script
Posh
9/2/2021
##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionThis is what your professor wants: “Can you please send a revised memo this week defining the empirical models you’ll run, and include a few baseline regressions (with regressions clearly defined), and include results for the simple models, including the correlation between violence and mop-up operations, with/without controls, and between lagged violence and mop-up operations, with/without controls.”
#Correlation between violence and mop-up
mop_viol_cor <- cor(conflict_data$REBEL.b, conflict_data$GOV_MOP.b) #a single correlation estimate
mop_viol_cor^2 #this is the proportion of mop-up operation variance shared with violence variance - R^2## [1] 0.08059437cor.test(conflict_data$REBEL.b, conflict_data$GOV_MOP.b) #with probability and other values##
## Pearson's product-moment correlation
##
## data: conflict_data$REBEL.b and conflict_data$GOV_MOP.b
## t = 189.89, df = 411364, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.2810795 0.2866986
## sample estimates:
## cor
## 0.2838915#Since it's a time series data we need to remove trend in time so the correlation between the variables
#is not inflated. For this we use REBEL.b.t1 variable (t-1) and GOV_MOP.b.t1. If there is a trend, the correlation might be
#smaller in the second case
mop_lag_viol_cor <- cor.test(conflict_data$REBEL.b.t1, conflict_data$GOV_MOP.b.t1)
mop_lag_viol_cor$estimate^2 #variance shared between mop-up and lagged violence## cor
## 0.08094853#All the rest of the analysis will use timelagged
#now we need to account for the variances shared by control variables. This is done by partial #correlation, using psych package.
partial.r(data=conflict_data, x=c("REBEL.b.t1","GOV_MOP.b.t1"), y=c("POP","ROAD_5KM","ELEVATION"))## partial correlations
## REBEL.b.t1 GOV_MOP.b.t1
## REBEL.b.t1 1.00 0.27
## GOV_MOP.b.t1 0.27 1.00# the predictive and the dependent variables are binary, so we need logistic regression. We use
# glm() function for logistic regression
logreg <- glm(REBEL.b.t1 ~ GOV_MOP.b.t1, data = conflict_data, family = binomial()) #without controls. #family specifies what kind of data we use for the regression
logreg_ctrl <- glm(REBEL.b.t1 ~ GOV_MOP.b.t1 + ELEVATION + POP + ROAD_5KM, data = conflict_data,
family = binomial()) #regression with controls
#let's see the summary of both models
summary(logreg)##
## Call:
## glm(formula = REBEL.b.t1 ~ GOV_MOP.b.t1, family = binomial(),
## data = conflict_data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.8850 -0.0788 -0.0788 -0.0788 3.3986
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -5.77215 0.02821 -204.58 <2e-16 ***
## GOV_MOP.b.t1 5.03688 0.06509 77.39 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 21912 on 407332 degrees of freedom
## Residual deviance: 18741 on 407331 degrees of freedom
## (4033 observations deleted due to missingness)
## AIC: 18745
##
## Number of Fisher Scoring iterations: 8summary(logreg_ctrl)##
## Call:
## glm(formula = REBEL.b.t1 ~ GOV_MOP.b.t1 + ELEVATION + POP + ROAD_5KM,
## family = binomial(), data = conflict_data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4204 -0.0985 -0.0675 -0.0449 4.0563
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.105e+00 1.305e-01 -46.794 < 2e-16 ***
## GOV_MOP.b.t1 4.263e+00 7.048e-02 60.494 < 2e-16 ***
## ELEVATION -1.079e-03 5.726e-05 -18.841 < 2e-16 ***
## POP 3.294e-04 1.046e-05 31.493 < 2e-16 ***
## ROAD_5KM 9.071e-01 1.296e-01 7.002 2.53e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 21907 on 406726 degrees of freedom
## Residual deviance: 17383 on 406722 degrees of freedom
## (4639 observations deleted due to missingness)
## AIC: 17393
##
## Number of Fisher Scoring iterations: 9#calculate chi statistic for model evaluation
logreg_chi <- logreg$null.deviance - logreg$deviance
#calculate df for chi
logregdf <- logreg$df.null - logreg$df.residual
#calculate probability for chi square
logreg_chisq_p <- 1 - pchisq(logreg_chi, logregdf)
#odds ratio for the model
exp(logreg$coefficients)## (Intercept) GOV_MOP.b.t1
## 3.113064e-03 1.539884e+02#calculate confidance intervals
exp(confint(logreg))## Waiting for profiling to be done...## 2.5 % 97.5 %
## (Intercept) 2.944072e-03 3.288417e-03
## GOV_MOP.b.t1 1.354689e+02 1.748523e+02#SECOND REGRESSION MODEL
#odds ratio for all the variables instead of log(odds ratio), for no reason
exp(logreg_ctrl$coefficients)## (Intercept) GOV_MOP.b.t1 ELEVATION POP ROAD_5KM
## 0.002231538 71.047974136 0.998921792 1.000329458 2.477110982#check the difference in deviance for the two models
diff_chi <- logreg$deviance - logreg_ctrl$deviance #the second model has a smaller deviance, #which is better
df_diff <- logreg$df.residual - logreg_ctrl$df.residual
diff_prob <- 1 - pchisq(diff_chi, df_diff)
#Let's check if the second model is the better one by calculating pseudo R^2.
#We subtract model deviance from the null deviance and divide by null deviance. Let's do this for the both models and then compare them. The model with the bigger R^2 is better as it #explains more variation in the dependent variable
#R^2 for the first model
m1_r2 <- (logreg$null.deviance - logreg$deviance)/logreg$null.deviance
#R^2 for the second model (with the controls)
m2_r2 <- (logreg_ctrl$null.deviance - logreg_ctrl$deviance)/logreg_ctrl$null.deviance
#The second model explains 20% of the variation while the first one explains the 14% of the #variation. For this reason the second model is a better model
#visualise logreg
conflict_data %>%
ggplot(aes(x = GOV_MOP.b.t1, y = REBEL.b.t1)) +
geom_jitter(alpha = 0.3, colour = "red3") +
theme_base() +
theme(plot.background = element_blank())## Warning: Removed 4033 rows containing missing values (geom_point).
model one - without controls Model deviance for this regression is less than baseline deviance, which means the model (regression) fits the data better than only using the intercept. Model x^2 (chi square) is 3170 with degrees of freedom 1. Including government operations as a predictor variable significantly improved the model (p = 0). This model shows that government intervention is a significant predictor of violent events, b = 5.04, z = 77.39, p < .001. Odds ratio is 153.9 (CI 135.5, 174.9), so if there is a government intervention it is highly likely that violence will occur, according to this model.
model two - with controls (population, accessibility and elevation) After adding three variables as controls in the regression model, b coefficient for government mop-up operation slightly decreased to 4.26. All the control variables significantly predict the violence but their coefficients are very low (elevation - -0.00, population - 0.00, accessibility - 0.91). Difference between the deviances of the models with and without control variables is statistically significant. Chi-square is 1372.05, df - 615, p = 0. The second model explains 20% of the variation while the first one explains the 14% of the variation. For this reason the second model is a better model. Odds ratio for the second model is a bit different. In this case government operations has an odds ratio of 71, almost half of that from the first model and accessible roads has an odds ratio of 2.5. This indicates, that including roads as a controlling variable might be beneficial for the theory.
##ლოგისტიკური რეგრესია
ლოგისტიკური რეგრესიის უკან ბევრი თეორიული კონცეპტი და ფორმულა დგას და ყველაფერი არც მე მესმის და შენც ვერ გაიგებ ეხლა დიდი დრო თუ არ დაუთმე. მოდელის შედეგები რომ გაიგო, რამდენიმე რაღაცას გავუსვათ ხაზი, რომ უფრო გასაგები იყოს.
პირველი, რახან შენი მთავარი დამოუკიდებელი და დამოკიდებული ცვლადი ბინარულია (0 და 1 პასუხებია მხოლოდ) შენი ჰიპოთეზა ერთის რაოდენობის ზრდასთან ერთად მეორეს შემთხვევების რაოდენობას ზრდას ვერ გაზომავს. მთავრობის ჩარევა რომ მატულობს რამდენი შემთხვევით მეტი გვაქვს ძალადობრივი ქმედება ამას ვერ გავიგებთ. ჩვენ გავიგებთ ერთი თვის განმავლობაში მთავრობა რომ ერევა, ერთხელ მაინც, რამდენად იზღდება იმის ალბათობა რომ ძალადობრივი ქმედება ერთხელ მაინც მოხდება. კი და არა პასუხები გვაქვს.
რეგრესიის წრფის გავლება ვერ მოხდება y ღერძე რახან მარტო ორი ნიშნული გვაქვს - 0 და 1. მაგიტომ გამოიყენება S ფორმის მრუდი. ეს მრუდი კი და არა პასუხების ნაცვლად უკვე გვეუბნება რამდენია იმის ალბათობა რომ ძალადობა მოხდეს. ნებისმიერი რიცხვი შეიძლება იყოს 0-დან ერთამდე. მაგრამ იმის გამო, რომ გვინდა რომ წრფივი რეგრესიის მსგავსად გამოვთვალოთ ლოგისტიკურის კოეფიციენტები, y ღერძი ისე უნდა შევცვალოთ, რომ წრფივი განტოლება მივიღოთ. ამისთვის ალბათობები გადაყავთ ლოგარითმულ სკალაზე და თან ალბათობების ფარდობას (ratio) იყენებენ პირდაპირ ალბათობის ნაცვლად. მაგალითად, ჩვენი მონაცემები რომ ავიღოთ და ერთი რომელიმე კონკრეტული მონაცემის ალბათობა 0.7 რომ იყოს, ლოგარითმულ (log odds) ღერძზე გადავა როგორც log(0.7/(1-0.7)), რაც 0.83-ია. ესე ვიღებთ სიმეტრიულ მონაცემებს. მარტივად რომ გითხრა, ლოგარითმი იმისთვის გამოიყენება, რომ როგორც ალბათობის სკალაზე 0-სა და 0.5-ს ან 0.5-სა და 1-ს შორის რიცხვები თანაბრადაა მოცემული, იგივე სიმეტრიულობით გვქონდეს რიცხვები წრფივი რეგრესიის სკალაზეც. ოღონდ 0-სა და 1-ს შორის მოქცეული რიცხვების ნაცვლად, რიცხვები მინუს უსასრულობიდან პლიუს უსასრულობამდე მიდის. ეს იმიტომ მოგიყევი, რომ ლოგისტიკურ ფლოთზე რომ ალბათობები წერია, ამის გამოთვლა სხვა y ღერძზე ხდება რეალურად (ლოგარითმულზე) და მერე გადმოგაქვს ალბათობის y ღერძზე. და კიდე უფრო მთავარი ისაა, რომ მოდელირებას რომ აკეთებ R-ში ან სხვა ნებისმიერ პროგრამაში, ლოგარითმულის პასუხებს გაძლევს და არა ალბათობებს. ალბათობებს მერე ითვლი. აქვე, ალბათობა რომ გამოთვალო, log(odds) კოეფიციენტი exp() ფორმულაში უნდა ჩასვა (it exponentiates it), მაგრამ აღარ არის ეს საჭირო.
მოდელის აღწერას რომ ნახულობ summary-ით, ესე უნდა გაიგო პასუხები (ბარემ ჩვენს მაგალითზე განვიხილოთ):
თავიდან გაჩვენებს უბრალოდ რა მოდელი გაუშვი (Call). ამას მოყვება deviance residuals აღწერა, რომელიც გეუბნება შენი რეალური მონაცემები რამდენადაა დაშორებული შენს მოდელს (დაახლოებით) და უკეთესი რომ აქ ყველა რიცხვი 0-თან ახლოს იყოს, ანუ ძალიან არ იყოს დაშ`რებული მოდელს.
ამის შემდეგ მოდის ყველაზე მნიშვნელოვანი ნაწილი - კოეფიციენტები. Intercept არის y ღერძზე წრფის გადაკვეთის ადგილი, ანუ x რომ 0-ს უდრის რას უდრის ამ დროს y. ჩვენ შემთხვევაში - რას უდრის ძალადობა, როცა მთავრობის ჩარევა 0-ის ტოლია. ჩვენ პირველ მოდელში ეს უდრის -5.77570-ს. st. error გეუბნება რეალობაში, ანუ მთავარ პოპულაციაში, საშუალოდ რამდენით განსხვავდება ეს რიცხვი პოპულაციის რიცხვისგან. ანუ ბევრი შერჩევა რომ აიღო და იმათში გაზომო ეს კოეფიციენტი საშუალოდ რა განსხვავება იქნება შენს ამ ერთ შერჩევაზე დაფუძნებით მიღებულ კოეფიციენტსა და სხვა ბევრ შერჩევაში მიღებულ კოეფიციენტებს შორის. Z სტატისტიკა გეუბნება 0-სგან განსხვავდება თუ - არა შენი კოეფიციენტი. და p გეუბნება ეს განსხვავება სტატისტიკურად მნიშვნელოვანია თუ - არა (95% ნდობის ინტერვალით). კოეფიციენტი 0-გან დიდად განსხვავებული რომ იყოს, p თუ მაღალია (>0.05) ნიშნავს რომ ეს შენი კოეფიციენტი ‘რეალურად’ შეიძლება არ განსხვავდებოდეს ნულისგან. და მთელი შენი იდეა რომ ერთი ცვლადი მეორეს პროგნოზირებს და გეუბნება რომ ამდენი და ამდენით იზრდება ან კლებულობს დამოკიდებული ცვლადი, რეალურად ეგრე არ არის. ჩვენ შემთხვევაში intercept-იც სტატისტიკურად მნიშვნელოვანია. ანუ როცა მთავრობა არ ერევა ძალადობის კოეფიციენტი მაგდენს უდრის.
მეორე კოეფიციენტი უკვე გეუბნება რამდენით შეიცვლება ძალადობა (log(odds) სკალაზე ვართ) თუკი ჩაურევლობიდან მთავრობის ჩარევაზე გადავალთ, ანუ 0-დან 1-ზე გადავალთ. ჩვენი პირველი მოდელი გვეუბნება, რომ ეს რიცხვი 5.03830-ია, ანუ ამდენით იზრდება. რახან ჩაურევლობის დროს ძალადობის მაჩვენებელი -5.77570-ია, ამას უნდა დავუმატოთ 5.03830 და მივიღებთ -0.7374-ს. დანარჩენი სტატისტიკები იგივეა რაც intercept-ზე.
ბოლო მნიშვნელოვანი რიცხვები მოდელის აღწერაში არის null.deviance და residual deviance. ამის უკან დგას მარტივი და მნიშვნელოვანი იდეა - თუ არსებობს იდეალური მოდელი, რომელიც სრული სიზუსტით პროგნოზირებს ჩვენს მონაცემებს, ჩვენ მიერ შექმნილი მოდელი რამდენად ახლოს შეიძლება იყოს მაგ იდეალურ მოდელთან. აი დაახლოებით ამას გეუბნება ეს რიცხვები. ამ იდეალურ მოდელს ადარებ ორ მოდელს - ერთს, რომელშიც ერთადერთი კოეფიციენტი არის საშუალო (null deviance), ანუ მარტო საშუალოს დახმარებით პროგნოზირებ მონაცემებს და მეორეს, სადაც უფრო მეტი პარამეტრი გაქვს (residual deviance). სადაც უფრო მეტი პარამეტრი შემოგაქვს (ჩვენ შემთხვევაში მთავრობის ჩარევის ფაქტორი) უფრო ნაკლებად უნდა იყოს დაშორებული იდეალურ მოდელს, ვიდრე ის მოდელი, სადაც მარტო საშუალო გაქვს. ანუ residual deviance რიცხვი უფრო ნაკლები უნდა იყოს, ვიდრე null deviance. ჩვენთან ესეა ორივე მოდელში, ამიტომ შეგვიძლია ვთქვათ, რომ საშუალოიან მოდელზე უკეთესი მოდელი შევქმენით.
ბოლო საკითხია ორი მოდელის შედარება. საკონტროლო ცვლადები შემოგვაქვს იმისთვის, რომ გავიგოთ სხვა ფაქტორებიც მოქმედებს თუ - არა დამოკიდებულ ცვლადზე. და რომ გავიგოთ ეს მოდელი უფრო კარგია თუ პირველი, სადაც საკონტროლო ცვლადებს არ ვიყენებთ. ერთი, რახან ვნახეთ, რომ მისასვლელობა სტატისტიკურად მნიშვნელოვნად პროგნოზირებს ძალადობას, ესეც მნიშვნელოვანია. მაგრამ R კვადრატს თუ გამოვთვლით ორივე მოდელისთვის (რაც გვეუბნება დამოკიდებულ ცვლადში ვარიაციის რა პროცენტს გვიხსნის მოდელი) გავიგებთ აუმჯობესებს თუ - არა და თან რამდენად, დამატებითი ცვლადის შემოტანა ჩვენს მოდელში. როგორც გამოვთვალეთ ზემოთ, ჩანს, რომ მისასვლელობის (გამოხატული გზის 5 კილომეტრიანი დაშორებით, ჩვენს მონაცემებში) დამატება მოდელში, აუმჯობესებს ჩვენს მოდელს - ეფექტის ზომას ზრდის, ვარიაციის მეტ პროცენტს გვიხსნის.
Thank you for sending, this is helpful to see.
There are more correlation tests in the memo than necessary for the final paper. With a simple observational correlation like this one, a few simple tests will be fine. This is true overall as well: the final document should include a few choice models, with more space devoted to careful discussion and interpretation.
I’d start with a simple: require(estimatr) lm_robust(violence ~ lagged_mop_up)
It was not clear to me in the memo if conflict_data$GOV_MOP.b is lagged MOP up operations. Note that you should lag your “treatment” (mop-up operations at t-1) and see how it correlates with violence in period t.
How many observations are there for each locale? If more than two time periods, you should be time-series regression with the plm() command. Please message me to let me know. I will be happy to send example code.
Then move on to lm_robust(violence ~ lagged_mop_up + controls)
In the document, please be clear about which control variables you include. This will be important for interpretation.
We’re just looking for the sign of the estimate and the statistical significance here, so the intuition from OLS is more useful, in my opinion, than the less intuitive glm logit/probit. The magnitude may also be of interest but it’s less important since the model is clearly not credible. This model simply shows the statistical relationship, and no evidence of the causal relationship.
The next steps will be to (1) the theory sections still needs discussion of COIN lit, like Biddle and Friedman. Some of the questions I asked above should be addressed in the theory section:
“Why do you predict that mop-up operations increase insurgent violence? What are the theoretical mechanisms linking mop-up operations and insurgent violence? What nodes lie between the two, and how do they aggregate such that COIN has a counterintuitive positive impact on violence? These mechanisms will be important for making theoretical predictions about cases when COIN may decrease violence, and when it may increase violence (such as your are predicting will be the case here). The existing literature, Biddle, Friedman, et al., will be the best place to turn to for offer theoretical statements about nodes between COIN and insurgent violence”
- Discuss the correlational results with a recognition of their limitations.
From earlier: “Think about potential confounds. The most obvious is reverse causation. If the government carries out COIN where there is more violence, we could observe a positive correlation between COIN and violence, even if COIN reduces violence. You can explore this threat to inference with placebo tests in the paper (see material from class and the textbook on placebo tests).”
- Then after that try matching on observable covariates: I see this as possible given your data, and a good way to use the methods from class for an extension
From earlier: “To provide more credible estimates of the relationship between mop-up operations and violence than naive regression, use methods from the class. Matching seems like a potentially useful approach for comparing areas with mop-up and without mop-up. A key part of the paper would be testing the suitability of the method, using the tests of design we ran in class.”
Which covariates do you have that you could match on?
In the class notes on matching please pay particular attention to the section on testing for the performance of matching. The important first step is to provide tests of whether the matching actually is a credible strategy or not.
Then you’ll want to move on to reestimating the relationship between lagged mop up operations and violence using the matched data. The discussion should include: 1) What is the sign, magnitude, and p value of relationship? 2) How wide is the 95% confidence interval? Is this a precise estimate or noisy estimate? 3) How does this estimate differ from the observational relationship? What does this imply about the potential confounds and reverse causation mentioned above? 4) What forms of confounding might still exist, even with matching?
Please keep me in the loop as you run these tests. Seen by Mariam Mgebrishvili at 9:34 AM _______________________________________________________________________________________
მივყვეთ ინსტრუქციას და ხელახლა გავუშვათ ახალი კოდი.
library(estimatr)
#regression without controls
rob_lm <- lm_robust(REBEL.b ~ GOV_MOP.b.t1, conflict_data)
#regression with control variables
rob_lm_ctrl <- lm_robust(REBEL.b ~ GOV_MOP.b.t1 + ELEVATION + POP + ROAD_5KM, data = conflict_data)
#visualise regression model (with controls)
plot_summs(rob_lm_ctrl, scale = TRUE)## Loading required namespace: broom.mixed
#fixed effect regression using plm
#creating plm data frame first
plm_conflict <- pdata.frame(conflict_data, index = c("NAME","YRMO"))## Warning in pdata.frame(conflict_data, index = c("NAME", "YRMO")): duplicate couples (id-time) in resulting pdata.frame
## to find out which, use, e.g., table(index(your_pdataframe), useNA = "ifany")#remove duplicates otherwise plm() won't work
plm_conflict$unique_id <- paste(plm_conflict$NAME, plm_conflict$YRMO)
plm_conflict$duplicate <- duplicated(plm_conflict$unique_id)
plm_conflict <- subset(plm_conflict, duplicate == "FALSE")
#calculate plm regression (take into account this takes some time)
model_plm <- plm(REBEL.b ~ GOV_MOP.b.t1 + diff(ELEVATION) + diff(POP) + diff(ROAD_5KM) + as.factor(NAME) + as.factor(YRMO), data = plm_conflict)
#plot the results for more clarity
plot_summs(model_plm)## Loading required namespace: broom.mixed
#let's see visually how treatment is spread across the dataset
DisplayTreatment(unit.id = "NAME", time.id = "YRMO",
xlab = "Year-Month", ylab = "Location",
treatment = "GOV_MOP.b", data = conflict_data,
hide.x.axis.label = TRUE, hide.y.axis.label = TRUE,
dense.plot = TRUE)
It looks the data is unbalanced with most of the locations and time there is no government intervention (red denotes to intervention at time t and location i [1], blue the absence of it [0]).
Before we do the estimation, we need to create matched sets across time.
#converting plm_datafram to regular data.frame
plm_conflict <- as.data.frame(plm_conflict)
#again, we are using pml_conflict because it doesn't contains duplicate rows
#changing class for the variables to use in the PanelMatch function
plm_conflict$CID <- as.integer(plm_conflict$CID)
plm_conflict$YRMO <- as.integer(plm_conflict$YRMO)
plm_conflict$GOV_MOP.b <- as.integer(plm_conflict$GOV_MOP.b)
plm_conflict$REBEL.b <- as.integer(plm_conflict$REBEL.b)
#creating matched sets
PM.results.none <- PanelMatch(lag = 4, time.id = "YRMO", unit.id = "CID",
treatment = "GOV_MOP.b", refinement.method = "none",
data = plm_conflict, match.missing = TRUE,
size.match = 5, qoi = "att", outcome.var = "REBEL.b",
lead = 0:4, forbid.treatment.reversal = FALSE,
use.diagonal.variance.matrix = TRUE)
#let's highlight a single match in the plot
#single out a match
mset <- PM.results.none$att[1]
# Use the DisplayTreatment function to visualize the
# treated unit and matched controls.
DisplayTreatment(unit.id = "CID",
time.id = "YRMO", legend.position = "none",
xlab = "year", ylab = "Location Code",
treatment = "GOV_MOP.b", data = plm_conflict,
matched.set = mset, # this way we highlight the particular set
show.set.only = TRUE)
We can ‘refine’ the matched data and use covariates for matching.
#note - needs to be specified whether we need to use lagged variables
PM.results.mah <- PanelMatch(lag = 4, time.id = "YRMO", unit.id = "CID",
treatment = "GOV_MOP.b", refinement.method = "mahalanobis",
data = plm_conflict, match.missing = TRUE,
covs.formula = ~ ELEVATION + POP + ROAD_5KM + ELECTIONS,
size.match = 5, qoi = "att", outcome.var = "REBEL.b",
lead = 0:4, forbid.treatment.reversal = FALSE,
use.diagonal.variance.matrix = TRUE)
#we can try a different method for refinement to see if the covariates are balanced.
#Extracting matched.set object from the result
msets.maha <- PM.results.mah$att
#see the results
print(msets.maha)## CID YRMO matched.set.size
## 1 22 33 4007
## 2 23 7 4003
## 3 23 14 3976
## 4 23 40 3981
## 5 23 44 8
## 6 23 47 11
## 7 23 51 4
## 8 23 54 11
## 9 23 57 13
## 10 23 75 3989
## 11 37 88 3978
## 12 43 85 3972
## 13 51 60 3946
## 14 61 76 3985
## 15 66 15 3971
## 16 67 52 3947
## 17 67 85 3972
## 18 76 57 3940
## 19 80 17 3964
## 20 166 10 4000
## 21 166 17 3964
## 22 166 40 3981
## 23 166 80 3974
## 24 166 96 4005
## 25 205 37 3986
## 26 206 7 4003
## 27 206 12 4000
## 28 206 25 3981
## 29 206 38 3983
## 30 206 46 3985
## 31 206 48 3
## 32 206 51 3
## 33 206 56 3940
## 34 206 75 3989
## 35 207 8 4002
## 36 207 20 3980
## 37 207 35 3997
## 38 207 43 3985
## 39 207 51 3954
## 40 207 58 3943
## 41 208 34 3998
## 42 208 60 3946
## 43 208 69 3988
## 44 209 7 4003
## 45 209 11 3
## 46 209 14 1
## 47 209 16 23
## 48 209 18 1
## 49 209 20 1
## 50 209 26 3981
## 51 209 32 4012
## 52 209 37 3986
## 53 209 40 11
## 54 209 50 3967
## 55 209 61 3951
## 56 209 65 2
## 57 209 72 3996
## 58 209 83 3977
## 59 212 62 3960
## 60 212 76 3985
## 61 212 82 3977
## 62 212 86 8
## 63 227 11 4000
## 64 227 14 1
## 65 227 21 3974
## 66 227 23 6
## 67 227 43 3985
## 68 227 58 3943
## 69 227 64 3984
## 70 227 78 3971
## 71 234 20 3980
## 72 335 8 4002
## 73 335 14 2
## 74 335 17 2
## 75 335 31 4013
## 76 335 38 3983
## 77 335 48 3989
## 78 335 52 1
## 79 335 54 11
## 80 335 56 1
## 81 335 63 4
## 82 335 67 2
## 83 335 75 3989
## 84 335 80 3974
## 85 374 84 3976
## 86 410 13 3993
## 87 410 25 3981
## 88 410 38 3983
## 89 410 51 3954
## 90 410 60 13
## 91 437 14 3976
## 92 437 17 2
## 93 437 23 3981
## 94 437 37 3986
## 95 437 39 14
## 96 437 48 3989
## 97 437 52 1
## 98 437 59 3946
## 99 437 62 11
## 100 437 66 7
## 101 437 68 8
## 102 437 75 3989
## 103 437 79 8
## 104 450 51 3954
## 105 450 83 3977
## 106 455 65 3991
## 107 467 18 3962
## 108 467 21 1
## 109 467 34 3998
## 110 467 50 3967
## 111 467 96 4005
## 112 480 56 3940
## 113 500 97 4008
## 114 503 57 3940
## 115 503 75 3989
## 116 522 79 3972
## 117 523 35 3997
## 118 523 52 3947
## 119 523 54 11
## 120 523 66 3988
## 121 523 68 8
## 122 553 87 3968
## 123 615 82 3977
## 124 615 86 8
## 125 640 14 3976
## 126 640 58 3943
## 127 640 76 3985
## 128 640 84 3976
## 129 657 25 3981
## 130 657 62 3960
## 131 668 17 3964
## 132 668 34 3998
## 133 676 68 3987
## 134 678 15 3971
## 135 678 25 3981
## 136 678 96 4005
## 137 699 50 3967
## 138 714 56 3940
## 139 714 58 9
## 140 714 63 3974
## 141 714 66 3
## 142 714 68 8
## 143 714 75 3989
## 144 714 84 3976
## 145 729 50 3967
## 146 772 49 3988
## 147 892 60 3946
## 148 892 66 7
## 149 896 55 3939
## 150 896 66 3988
## 151 907 62 3960
## 152 907 64 9
## 153 962 44 3983
## 154 962 56 3940
## 155 962 65 3991
## 156 976 36 3996
## 157 976 53 3933
## 158 976 75 3989
## 159 976 77 10
## 160 976 80 12
## 161 1047 11 4000
## 162 1047 17 3964
## 163 1047 50 3967
## 164 1047 57 3940
## 165 1047 59 13
## 166 1047 82 3977
## 167 1049 92 4005
## 168 1052 17 3964
## 169 1054 21 3974
## 170 1054 44 3983
## 171 1054 49 3988
## 172 1054 52 10
## 173 1054 58 3943
## 174 1054 96 4005
## 175 1094 6 4000
## 176 1094 14 3976
## 177 1105 17 3964
## 178 1118 57 3940
## 179 1120 17 3964
## 180 1120 24 9
## 181 1120 39 3985
## 182 1120 51 3954
## 183 1123 53 3933
## 184 1139 53 3933
## 185 1154 50 3967
## 186 1180 10 4000
## 187 1180 13 4
## 188 1180 19 3982
## 189 1180 50 3967
## 190 1180 52 19
## 191 1180 55 12
## 192 1180 67 3992
## 193 1193 37 3986
## 194 1193 40 11
## 195 1193 47 3985
## 196 1193 50 5
## 197 1193 52 19
## 198 1193 54 3
## 199 1193 57 13
## 200 1193 81 3979
## 201 1206 31 4013
## 202 1214 13 3993
## 203 1214 15 8
## 204 1214 52 3947
## 205 1214 56 2
## 206 1214 82 3977
## 207 1245 64 3984
## 208 1279 7 5
## 209 1279 12 4000
## 210 1279 14 3
## 211 1279 18 27
## 212 1279 21 1
## 213 1279 23 1
## 214 1279 47 3985
## 215 1279 54 3933
## 216 1289 9 4000
## 217 1289 52 3947
## 218 1289 75 3989
## 219 1291 13 3993
## 220 1303 73 4004
## 221 1303 77 3
## 222 1335 47 3985
## 223 1335 50 5
## 224 1335 55 3939
## 225 1335 77 3975
## 226 1335 88 3978
## 227 1346 27 3982
## 228 1367 68 3987
## 229 1367 92 4005
## 230 1367 97 4008
## 231 1382 14 3976
## 232 1382 19 3982
## 233 1455 59 3946
## 234 1485 15 3971
## 235 1485 38 3983
## 236 1485 42 3
## 237 1485 49 3988
## 238 1502 13 3993
## 239 1502 15 8
## 240 1502 20 3980
## 241 1502 32 4012
## 242 1502 34 4
## 243 1502 53 3933
## 244 1502 57 12
## 245 1507 75 3989
## 246 1509 54 3933
## 247 1547 25 3981
## 248 1547 42 3995
## 249 1547 48 3989
## 250 1547 62 3960
## 251 1560 34 3998
## 252 1567 17 3964
## 253 1567 52 3947
## 254 1567 74 4001
## 255 1574 79 3972
## 256 1615 17 3964
## 257 1615 23 3981
## 258 1615 53 3933
## 259 1729 43 3985
## 260 1731 55 3939
## 261 1735 34 3998
## 262 1735 36 10
## 263 1735 56 3940
## 264 1735 85 3972
## 265 1735 94 4012
## 266 1747 57 3940
## 267 1756 21 3974
## 268 1756 47 3985
## 269 1756 52 3947
## 270 1756 60 3946
## 271 1767 42 3995
## 272 1793 87 3968
## 273 1802 8 4002
## 274 1802 14 3976
## 275 1802 16 23
## 276 1802 18 1
## 277 1802 26 3981
## 278 1802 31 4
## 279 1802 45 3989
## 280 1802 54 3933
## 281 1802 62 3960
## 282 1802 72 3996
## 283 1803 18 3962
## 284 1803 26 3981
## 285 1803 44 3983
## 286 1803 53 3933
## 287 1803 81 3979
## 288 1808 14 3976
## 289 1808 23 3981
## 290 1808 25 9
## 291 1808 31 4013
## 292 1808 42 3995
## 293 1808 45 8
## 294 1808 48 2
## 295 1808 52 1
## 296 1808 55 12
## 297 1826 61 3951
## 298 1852 13 3993
## 299 1852 25 3981
## 300 1856 44 3983
## 301 1856 50 3967
## 302 1864 5 3999
## 303 1876 14 3976
## 304 1876 20 3980
## 305 1876 58 3943
## 306 1876 83 3977
## 307 1880 31 4013
## 308 1883 82 3977
## 309 1883 90 3989
## 310 1926 20 3980
## 311 1926 27 3982
## 312 1926 37 3986
## 313 1926 48 3989
## 314 1926 62 3960
## 315 1926 65 11
## 316 1926 72 3996
## 317 1926 76 2
## 318 1926 83 3977
## 319 1927 7 4003
## 320 1927 14 3976
## 321 1927 19 3982
## 322 1927 66 3988
## 323 1927 92 4005
## 324 1929 14 3976
## 325 1929 23 3981
## 326 1929 39 3985
## 327 1929 55 3939
## 328 1960 87 3968
## 329 1968 50 3967
## 330 1968 66 3988
## 331 1970 18 3962
## 332 1973 37 3986
## 333 2010 18 3962
## 334 2010 86 3975
## 335 2010 88 4
## 336 2037 50 3967
## 337 2037 52 19
## 338 2065 42 1
## 339 2065 96 4005
## 340 2080 82 3977
## 341 2080 88 3978
## 342 2087 6 4000
## 343 2087 8 4
## 344 2087 10 1
## 345 2087 13 4
## 346 2087 21 1
## 347 2087 31 4013
## 348 2087 35 5
## 349 2087 38 3
## 350 2087 43 11
## 351 2087 49 3988
## 352 2087 56 2
## 353 2087 59 8
## 354 2087 62 4
## 355 2087 66 7
## 356 2087 68 8
## 357 2087 73 7
## 358 2087 75 3
## 359 2087 81 11
## 360 2087 83 7
## 361 2087 88 1
## 362 2087 94 4012
## 363 2128 43 3985
## 364 2175 53 3933
## 365 2181 8 4002
## 366 2181 14 3976
## 367 2181 20 3980
## 368 2181 49 3988
## 369 2181 52 10
## 370 2181 82 3977
## 371 2192 57 3940
## 372 2196 21 3974
## 373 2196 34 3998
## 374 2196 42 3995
## 375 2196 48 3989
## 376 2196 51 5
## 377 2196 54 11
## 378 2196 56 15
## 379 2196 60 6
## 380 2196 74 4001
## 381 2196 79 8
## 382 2196 81 7
## 383 2196 83 1
## 384 2196 88 7
## 385 2196 96 4005
## 386 2200 14 3976
## 387 2200 37 3986
## 388 2266 23 3981
## 389 2266 57 3940
## 390 2266 82 3977
## 391 2300 88 3978
## 392 2322 46 3985
## 393 2327 22 3979
## 394 2380 51 3954
## 395 2428 23 3981
## 396 2428 39 3985
## 397 2428 58 3943
## 398 2430 55 3939
## 399 2430 59 11
## 400 2430 77 3975
## 401 2430 84 3976
## 402 2437 9 4000
## 403 2437 12 4
## 404 2437 18 3962
## 405 2437 87 3968
## 406 2437 96 4005
## 407 2440 92 4005
## 408 2445 15 3971
## 409 2445 53 3933
## 410 2445 66 3988
## 411 2450 59 3946
## 412 2450 73 4004
## 413 2492 24 3985
## 414 2492 34 3998
## 415 2492 49 3988
## 416 2492 59 3946
## 417 2500 87 3968
## 418 2501 86 3975
## 419 2503 43 3985
## 420 2503 53 3933
## 421 2564 37 3986
## 422 2564 48 3989
## 423 2564 51 5
## 424 2568 65 3991
## 425 2568 85 3972
## 426 2578 30 4014
## 427 2578 35 3997
## 428 2578 51 3954
## 429 2585 79 3972
## 430 2588 54 3933
## 431 2588 68 3987
## 432 2589 50 3967
## 433 2611 80 3974
## 434 2615 15 3971
## 435 2615 36 3996
## 436 2615 68 3987
## 437 2615 75 3989
## 438 2615 89 3984
## 439 2709 10 4000
## 440 2714 57 3940
## 441 2745 48 3989
## 442 2745 51 3
## 443 2745 53 15
## 444 2767 81 3979
## 445 2784 22 3979
## 446 2784 55 3939
## 447 2784 68 3987
## 448 2848 40 3981
## 449 2848 55 3939
## 450 2863 34 3998
## 451 2863 37 9
## 452 2863 49 3988
## 453 2863 59 3946
## 454 2863 63 16
## 455 2863 67 2
## 456 2863 76 3985
## 457 2863 78 9
## 458 2863 84 3976
## 459 2863 88 9
## 460 2863 96 4005
## 461 2907 78 3971
## 462 2954 83 3977
## 463 2972 36 3996
## 464 2975 74 4001
## 465 2982 48 3989
## 466 2982 56 3940
## 467 2982 58 9
## 468 2982 70 3991
## 469 2982 73 3
## 470 2982 84 3976
## 471 2990 52 3947
## 472 2990 66 3988
## 473 2990 85 3972
## 474 2990 88 7
## 475 3044 67 3992
## 476 3048 77 3975
## 477 3056 49 3988
## 478 3056 79 3972
## 479 3282 79 3972
## 480 3285 78 3971
## 481 3315 75 3989
## 482 3315 79 8
## 483 3315 87 3968
## 484 3322 50 3967
## 485 3322 60 3946
## 486 3369 60 3946
## 487 3402 18 3962
## 488 3424 50 3967
## 489 3424 52 19
## 490 3424 74 4001
## 491 3424 77 8
## 492 3424 87 3968
## 493 3424 89 14
## 494 3424 94 4012
## 495 3424 96 3
## 496 3439 20 3980
## 497 3439 23 7
## 498 3439 27 10
## 499 3439 37 3986
## 500 3439 47 3985
## 501 3439 49 7
## 502 3439 52 10
## 503 3439 64 3984
## 504 3439 69 3988
## 505 3439 76 3985
## 506 3439 83 3977
## 507 3498 21 3974
## 508 3498 43 3985
## 509 3498 53 3933
## 510 3498 64 3984
## 511 3667 5 3999
## 512 3667 38 3983
## 513 3866 14 3976
## 514 3866 50 3967
## 515 3866 58 3943
## 516 3924 12 4000
## 517 3924 21 8
## 518 3924 24 6
## 519 3924 30 4014
## 520 3924 42 3995
## 521 3924 44 8
## 522 3924 49 2
## 523 3924 57 3940
## 524 3924 70 3991
## 525 3924 78 3971
## 526 3976 88 3978
## 527 4078 13 3993
## 528 4078 16 7
## 529 4078 22 3979
## 530 4078 48 3989
## 531 4078 53 3933
## 532 4078 59 3946
## 533 4078 70 3991
## 534 4088 89 3984
## 535 4126 36 3996
## 536 4137 59 3946
## 537 4137 75 3989
## 538 4137 93 4012
## 539 4285 83 3977
## 540 4287 76 3985
## 541 4288 55 3939
## 542 4294 48 3989
## 543 4294 53 3933
## 544 4294 55 12
## 545 4294 58 12
## 546 4294 70 3991
## 547 4294 76 3985
## 548 4294 79 8
## 549 4294 81 7
## 550 4294 85 6
## 551 4294 88 7
## 552 4295 60 3946
## 553 4309 5 3999
## 554 4309 25 3981
## 555 4309 32 4012
## 556 4309 36 3
## 557 4309 59 3946
## 558 4310 74 4001
## 559 4313 75 3989
## 560 4315 76 3985
## 561 4358 44 3983
## 562 4358 49 2
## 563 4358 75 3989
## 564 4358 80 11
## 565 4358 83 5
## 566 4358 85 13
## 567 4358 94 4012
## 568 4358 96 3
## 569 4452 87 3968
## 570 4452 93 4012
## 571 4457 9 4000
## 572 4457 21 3974
## 573 4457 34 3998
## 574 4457 55 3939
## 575 4457 57 12
## 576 4457 69 3988
## 577 4457 80 3974
## 578 4457 93 4012
## 579 4458 46 3985
## 580 4458 51 3954
## 581 4458 59 11
## 582 4463 69 3988
## 583 4540 69 3988
## 584 4543 6 4
## 585 4543 14 3976
## 586 4543 19 3982
## 587 4543 23 5
## 588 4543 57 3940
## 589 4543 74 4001
## 590 4543 84 3976
## 591 4544 50 3967
## 592 4544 59 3946
## 593 4559 6 4000
## 594 4559 14 3976
## 595 4559 43 3985
## 596 4559 50 3967
## 597 4559 61 3951
## 598 4559 81 3979
## 599 4731 39 3985
## 600 4731 59 3946
## 601 4731 79 3972
## 602 4734 82 3977
## 603 4760 78 3971
## 604 4760 87 3968
## 605 4800 20 3980
## 606 4800 23 7
## 607 4800 32 4012
## 608 4800 50 3967
## 609 4800 53 18
## 610 4800 55 12
## 611 4800 58 12
## 612 4800 66 3988
## 613 4800 69 7
## 614 4803 8 4002
## 615 4808 42 3995
## 616 4817 64 3984
## 617 4817 72 1
## 618 4817 81 3
## 619 4817 84 7
## 620 4817 87 10
## 621 4817 89 14
## 622 4817 92 4
## 623 4843 50 3967
## 624 4853 21 3974
## 625 4853 26 3981
## 626 4853 42 3995
## 627 4853 48 3989
## 628 4853 53 2
## 629 4853 55 12
## 630 4853 66 7
## 631 4853 68 8
## 632 4853 75 3989
## 633 4853 80 1
## 634 4853 93 5
## 635 4853 98 3
## 636 4857 70 3991
## 637 4857 80 3974
## 638 4889 56 3940
## 639 4889 60 9
## 640 4889 68 3987
## 641 4969 33 4007
## 642 4969 36 2
## 643 4969 53 3933
## 644 4969 87 3968
## 645 5054 83 3977
## 646 5158 9 4000
## 647 5158 43 3985
## 648 5158 54 3933
## 649 5158 59 3946
## 650 5158 78 3971
## 651 5200 51 3954
## 652 5200 54 11
## 653 5200 56 15
## 654 5306 43 3985
## 655 5338 10 1
## 656 5338 12 1
## 657 5338 24 1
## 658 5338 30 4014
## 659 5338 38 3983
## 660 5338 44 3983
## 661 5338 49 3988
## 662 5338 59 3946
## 663 5338 64 14
## 664 5338 74 4001
## 665 5360 87 3968
## 666 5426 22 3979
## 667 5426 52 3947
## 668 5426 56 12
## 669 5431 8 4002
## 670 5431 12 5
## 671 5431 80 3974
## 672 5467 43 3985
## 673 5467 53 3933
## 674 5472 81 3979
## 675 5495 9 4000
## 676 5495 12 4
## 677 5495 19 3982
## 678 5495 23 1
## 679 5495 32 4012
## 680 5495 36 3
## 681 5495 39 10
## 682 5495 43 11
## 683 5495 50 3967
## 684 5495 52 19
## 685 5495 56 2
## 686 5495 59 8
## 687 5495 61 11
## 688 5495 65 4
## 689 5495 67 8
## 690 5495 74 4001
## 691 5495 78 7
## 692 5495 86 3975
## 693 5524 25 3981
## 694 5524 58 3943
## 695 5524 83 3977
## 696 5531 76 3985
## 697 5538 14 3976
## 698 5538 23 3981
## 699 5700 19 3982
## 700 5700 23 5
## 701 5700 46 3985
## 702 5700 53 3933
## 703 5700 57 12
## 704 5702 13 3993
## 705 5714 13 3993
## 706 5714 18 3962
## 707 5714 62 3960
## 708 5763 71 3995
## 709 5844 82 3977
## 710 5899 53 3933
## 711 5947 6 4000
## 712 5947 9 4
## 713 5947 13 1
## 714 5947 17 4
## 715 5947 19 9
## 716 5947 23 3
## 717 5947 28 3991
## 718 5947 40 3981
## 719 5947 43 8
## 720 5947 50 3967
## 721 5947 52 19
## 722 5947 74 4001
## 723 5947 78 7
## 724 5947 96 4005
## 725 6014 9 4000
## 726 6014 14 3976
## 727 6014 20 3980
## 728 6014 26 3981
## 729 6014 36 3996
## 730 6014 43 3985
## 731 6014 51 3954
## 732 6014 55 10
## 733 6014 58 12
## 734 6067 54 3933
## 735 6133 58 3943
## 736 6133 96 4005
## 737 6159 51 3954
## 738 6212 64 3984
## 739 6238 25 3981
## 740 6238 61 3951
## 741 6285 37 3986
## 742 6327 13 3993
## 743 6327 19 3982
## 744 6327 34 3998
## 745 6327 40 3981
## 746 6327 48 3989
## 747 6327 51 3
## 748 6327 54 11
## 749 6327 62 3960
## 750 6327 71 3995
## 751 6327 96 4005
## 752 6354 25 3981
## 753 6354 54 3933
## 754 6354 92 4005
## 755 6444 97 4008
## 756 6447 71 3995
## 757 6462 30 4014
## 758 6462 36 3996
## 759 6462 39 10
## 760 6462 49 3988
## 761 6462 54 3933
## 762 6462 57 13
## 763 6462 63 16
## 764 6462 69 3988
## 765 6462 74 4001
## 766 6462 85 3972
## 767 6462 88 7
## 768 6506 75 3989
## 769 6506 81 3979
## 770 6514 77 3975
## 771 6537 5 3999
## 772 6537 8 4
## 773 6537 43 3985
## 774 6537 84 3976
## 775 6582 13 3993
## 776 6582 39 3985
## 777 6582 52 3947
## 778 6582 57 5
## 779 6582 68 3987
## 780 6582 70 6
## 781 6583 51 3954
## 782 6588 15 3971
## 783 6588 18 7
## 784 6588 21 8
## 785 6588 23 6
## 786 6588 25 3
## 787 6588 39 3985
## 788 6588 44 3983
## 789 6588 48 2
## 790 6588 52 1
## 791 6588 62 3960
## 792 6588 85 3972
## 793 6640 13 3993
## 794 6640 26 3981
## 795 6640 44 3983
## 796 6640 51 3954
## 797 6640 54 11
## 798 6640 56 15
## 799 6640 77 3975
## 800 6641 6 4000
## 801 6641 8 4
## 802 6641 12 5
## 803 6641 16 1
## 804 6641 20 1
## 805 6641 23 7
## 806 6641 25 9
## 807 6641 40 3981
## 808 6641 42 9
## 809 6641 47 3985
## 810 6641 53 3933
## 811 6641 59 2
## 812 6666 73 4004
## 813 6667 14 3976
## 814 6667 52 3947
## 815 6670 17 3964
## 816 6670 26 3981
## 817 6670 60 3946
## 818 6670 77 3975
## 819 6687 62 3960
## 820 6687 66 7
## 821 6687 77 3975
## 822 6687 83 3977
## 823 6743 5 3999
## 824 6802 54 3933
## 825 6819 56 3940
## 826 6840 57 3940
## 827 6937 84 3976
## 828 7024 33 4007
## 829 7024 36 2
## 830 7024 87 3968
## 831 7035 18 3962
## 832 7046 10 4000
## 833 7046 14 2
## 834 7046 33 4007
## 835 7046 39 3985
## 836 7046 50 3967
## 837 7046 54 18
## 838 7046 58 16
## 839 7046 77 3975
## 840 7095 14 3976
## 841 7095 20 3980
## 842 7096 44 3983
## 843 7096 51 3954
## 844 7096 65 3991
## 845 7096 77 3975
## 846 7097 44 3983
## 847 7097 51 3954
## 848 7097 65 3991
## 849 7097 77 3975
## 850 7125 83 3977
## 851 7202 87 3968
## 852 7221 10 4000
## 853 7221 14 2
## 854 7221 33 4007
## 855 7221 39 3985
## 856 7221 50 3967
## 857 7221 54 18
## 858 7221 58 16
## 859 7221 77 3975
## 860 7226 44 3983
## 861 7236 40 3981
## 862 7236 57 3940
## 863 7236 82 3977
## 864 7236 84 9
## 865 7236 89 3984
## 866 7236 93 5
## 867 7249 14 3976
## 868 7314 85 3972
## 869 7345 18 3962
## 870 7345 25 3981
## 871 7345 43 3985
## 872 7345 56 3940
## 873 7367 54 3933
## 874 7401 37 3986
## 875 7404 45 3989
## 876 7449 83 3977
## 877 7568 44 3983
## 878 7568 67 3992
## 879 7583 50 3967
## 880 7601 60 3946
## 881 7601 90 3989
## 882 7798 65 3991
## 883 7798 85 3972
## 884 7860 78 3971
## 885 7876 47 3985
## 886 7897 87 3968
## 887 7902 83 3977
## 888 7946 6 7
## 889 7946 11 4000
## 890 7946 19 8
## 891 7946 21 5
## 892 7946 23 1
## 893 7946 27 10
## 894 7946 29 3
## 895 7946 34 1
## 896 7946 37 9
## 897 7946 47 3985
## 898 7946 49 7
## 899 7946 52 10
## 900 7946 54 11
## 901 7946 72 3996
## 902 7946 85 3972
## 903 7948 90 3989
## 904 8000 14 3976
## 905 8000 29 3996
## 906 8000 42 3995
## 907 8034 11 4000
## 908 8034 14 1
## 909 8034 25 3981
## 910 8034 37 3986
## 911 8034 40 11
## 912 8034 49 3988
## 913 8034 59 3946
## 914 8034 62 4
## 915 8034 71 3995
## 916 8035 35 3997
## 917 8035 60 3946
## 918 8057 33 4007
## 919 8057 36 2
## 920 8057 87 3968
## 921 8060 15 3971
## 922 8060 55 3939
## 923 8257 14 3976
## 924 8288 44 3983
## 925 8288 82 3977
## 926 8288 84 9
## 927 8422 51 3954
## 928 8422 60 3946
## 929 8422 66 7
## 930 8422 68 8
## 931 8422 76 3985
## 932 8422 94 4012
## 933 8432 82 3977
## 934 8502 37 3986
## 935 8512 65 3991
## 936 8513 52 3947
## 937 8513 77 3975
## 938 8519 25 3981
## 939 8519 51 3954
## 940 8519 54 11
## 941 8519 56 15
## 942 8599 17 3964
## 943 8599 19 9
## 944 8599 49 3988
## 945 8599 55 3939
## 946 8637 54 3933
## 947 8694 13 3993
## 948 8694 55 3939
## 949 335 25 0
## 950 410 57 0
## 951 1120 55 0
## 952 1214 19 0
## 953 2065 25 0
## 954 2065 28 0
## 955 2065 30 0
## 956 2065 74 0
## 957 2065 83 0
## 958 2065 86 0
## 959 2065 88 0
## 960 2087 17 0
## 961 2087 25 0
## 962 2087 77 0
## 963 2181 25 0
## 964 2196 85 0
## 965 3924 15 0
## 966 3924 18 0
## 967 4457 25 0
## 968 4458 55 0
## 969 4817 68 0
## 970 4853 60 0
## 971 4853 62 0
## 972 4853 82 0
## 973 5338 7 0
## 974 5338 14 0
## 975 5338 18 0
## 976 5431 84 0
## 977 5947 57 0
## 978 7946 14 0
## 979 8422 79 0summary(msets.maha)## $overview
## CID YRMO matched.set.size
## 1 22 33 4007
## 2 23 7 4003
## 3 23 14 3976
## 4 23 40 3981
## 5 23 44 8
## 6 23 47 11
## 7 23 51 4
## 8 23 54 11
## 9 23 57 13
## 10 23 75 3989
## 11 37 88 3978
## 12 43 85 3972
## 13 51 60 3946
## 14 61 76 3985
## 15 66 15 3971
## 16 67 52 3947
## 17 67 85 3972
## 18 76 57 3940
## 19 80 17 3964
## 20 166 10 4000
## 21 166 17 3964
## 22 166 40 3981
## 23 166 80 3974
## 24 166 96 4005
## 25 205 37 3986
## 26 206 7 4003
## 27 206 12 4000
## 28 206 25 3981
## 29 206 38 3983
## 30 206 46 3985
## 31 206 48 3
## 32 206 51 3
## 33 206 56 3940
## 34 206 75 3989
## 35 207 8 4002
## 36 207 20 3980
## 37 207 35 3997
## 38 207 43 3985
## 39 207 51 3954
## 40 207 58 3943
## 41 208 34 3998
## 42 208 60 3946
## 43 208 69 3988
## 44 209 7 4003
## 45 209 11 3
## 46 209 14 1
## 47 209 16 23
## 48 209 18 1
## 49 209 20 1
## 50 209 26 3981
## 51 209 32 4012
## 52 209 37 3986
## 53 209 40 11
## 54 209 50 3967
## 55 209 61 3951
## 56 209 65 2
## 57 209 72 3996
## 58 209 83 3977
## 59 212 62 3960
## 60 212 76 3985
## 61 212 82 3977
## 62 212 86 8
## 63 227 11 4000
## 64 227 14 1
## 65 227 21 3974
## 66 227 23 6
## 67 227 43 3985
## 68 227 58 3943
## 69 227 64 3984
## 70 227 78 3971
## 71 234 20 3980
## 72 335 8 4002
## 73 335 14 2
## 74 335 17 2
## 75 335 31 4013
## 76 335 38 3983
## 77 335 48 3989
## 78 335 52 1
## 79 335 54 11
## 80 335 56 1
## 81 335 63 4
## 82 335 67 2
## 83 335 75 3989
## 84 335 80 3974
## 85 374 84 3976
## 86 410 13 3993
## 87 410 25 3981
## 88 410 38 3983
## 89 410 51 3954
## 90 410 60 13
## 91 437 14 3976
## 92 437 17 2
## 93 437 23 3981
## 94 437 37 3986
## 95 437 39 14
## 96 437 48 3989
## 97 437 52 1
## 98 437 59 3946
## 99 437 62 11
## 100 437 66 7
## 101 437 68 8
## 102 437 75 3989
## 103 437 79 8
## 104 450 51 3954
## 105 450 83 3977
## 106 455 65 3991
## 107 467 18 3962
## 108 467 21 1
## 109 467 34 3998
## 110 467 50 3967
## 111 467 96 4005
## 112 480 56 3940
## 113 500 97 4008
## 114 503 57 3940
## 115 503 75 3989
## 116 522 79 3972
## 117 523 35 3997
## 118 523 52 3947
## 119 523 54 11
## 120 523 66 3988
## 121 523 68 8
## 122 553 87 3968
## 123 615 82 3977
## 124 615 86 8
## 125 640 14 3976
## 126 640 58 3943
## 127 640 76 3985
## 128 640 84 3976
## 129 657 25 3981
## 130 657 62 3960
## 131 668 17 3964
## 132 668 34 3998
## 133 676 68 3987
## 134 678 15 3971
## 135 678 25 3981
## 136 678 96 4005
## 137 699 50 3967
## 138 714 56 3940
## 139 714 58 9
## 140 714 63 3974
## 141 714 66 3
## 142 714 68 8
## 143 714 75 3989
## 144 714 84 3976
## 145 729 50 3967
## 146 772 49 3988
## 147 892 60 3946
## 148 892 66 7
## 149 896 55 3939
## 150 896 66 3988
## 151 907 62 3960
## 152 907 64 9
## 153 962 44 3983
## 154 962 56 3940
## 155 962 65 3991
## 156 976 36 3996
## 157 976 53 3933
## 158 976 75 3989
## 159 976 77 10
## 160 976 80 12
## 161 1047 11 4000
## 162 1047 17 3964
## 163 1047 50 3967
## 164 1047 57 3940
## 165 1047 59 13
## 166 1047 82 3977
## 167 1049 92 4005
## 168 1052 17 3964
## 169 1054 21 3974
## 170 1054 44 3983
## 171 1054 49 3988
## 172 1054 52 10
## 173 1054 58 3943
## 174 1054 96 4005
## 175 1094 6 4000
## 176 1094 14 3976
## 177 1105 17 3964
## 178 1118 57 3940
## 179 1120 17 3964
## 180 1120 24 9
## 181 1120 39 3985
## 182 1120 51 3954
## 183 1123 53 3933
## 184 1139 53 3933
## 185 1154 50 3967
## 186 1180 10 4000
## 187 1180 13 4
## 188 1180 19 3982
## 189 1180 50 3967
## 190 1180 52 19
## 191 1180 55 12
## 192 1180 67 3992
## 193 1193 37 3986
## 194 1193 40 11
## 195 1193 47 3985
## 196 1193 50 5
## 197 1193 52 19
## 198 1193 54 3
## 199 1193 57 13
## 200 1193 81 3979
## 201 1206 31 4013
## 202 1214 13 3993
## 203 1214 15 8
## 204 1214 52 3947
## 205 1214 56 2
## 206 1214 82 3977
## 207 1245 64 3984
## 208 1279 7 5
## 209 1279 12 4000
## 210 1279 14 3
## 211 1279 18 27
## 212 1279 21 1
## 213 1279 23 1
## 214 1279 47 3985
## 215 1279 54 3933
## 216 1289 9 4000
## 217 1289 52 3947
## 218 1289 75 3989
## 219 1291 13 3993
## 220 1303 73 4004
## 221 1303 77 3
## 222 1335 47 3985
## 223 1335 50 5
## 224 1335 55 3939
## 225 1335 77 3975
## 226 1335 88 3978
## 227 1346 27 3982
## 228 1367 68 3987
## 229 1367 92 4005
## 230 1367 97 4008
## 231 1382 14 3976
## 232 1382 19 3982
## 233 1455 59 3946
## 234 1485 15 3971
## 235 1485 38 3983
## 236 1485 42 3
## 237 1485 49 3988
## 238 1502 13 3993
## 239 1502 15 8
## 240 1502 20 3980
## 241 1502 32 4012
## 242 1502 34 4
## 243 1502 53 3933
## 244 1502 57 12
## 245 1507 75 3989
## 246 1509 54 3933
## 247 1547 25 3981
## 248 1547 42 3995
## 249 1547 48 3989
## 250 1547 62 3960
## 251 1560 34 3998
## 252 1567 17 3964
## 253 1567 52 3947
## 254 1567 74 4001
## 255 1574 79 3972
## 256 1615 17 3964
## 257 1615 23 3981
## 258 1615 53 3933
## 259 1729 43 3985
## 260 1731 55 3939
## 261 1735 34 3998
## 262 1735 36 10
## 263 1735 56 3940
## 264 1735 85 3972
## 265 1735 94 4012
## 266 1747 57 3940
## 267 1756 21 3974
## 268 1756 47 3985
## 269 1756 52 3947
## 270 1756 60 3946
## 271 1767 42 3995
## 272 1793 87 3968
## 273 1802 8 4002
## 274 1802 14 3976
## 275 1802 16 23
## 276 1802 18 1
## 277 1802 26 3981
## 278 1802 31 4
## 279 1802 45 3989
## 280 1802 54 3933
## 281 1802 62 3960
## 282 1802 72 3996
## 283 1803 18 3962
## 284 1803 26 3981
## 285 1803 44 3983
## 286 1803 53 3933
## 287 1803 81 3979
## 288 1808 14 3976
## 289 1808 23 3981
## 290 1808 25 9
## 291 1808 31 4013
## 292 1808 42 3995
## 293 1808 45 8
## 294 1808 48 2
## 295 1808 52 1
## 296 1808 55 12
## 297 1826 61 3951
## 298 1852 13 3993
## 299 1852 25 3981
## 300 1856 44 3983
## 301 1856 50 3967
## 302 1864 5 3999
## 303 1876 14 3976
## 304 1876 20 3980
## 305 1876 58 3943
## 306 1876 83 3977
## 307 1880 31 4013
## 308 1883 82 3977
## 309 1883 90 3989
## 310 1926 20 3980
## 311 1926 27 3982
## 312 1926 37 3986
## 313 1926 48 3989
## 314 1926 62 3960
## 315 1926 65 11
## 316 1926 72 3996
## 317 1926 76 2
## 318 1926 83 3977
## 319 1927 7 4003
## 320 1927 14 3976
## 321 1927 19 3982
## 322 1927 66 3988
## 323 1927 92 4005
## 324 1929 14 3976
## 325 1929 23 3981
## 326 1929 39 3985
## 327 1929 55 3939
## 328 1960 87 3968
## 329 1968 50 3967
## 330 1968 66 3988
## 331 1970 18 3962
## 332 1973 37 3986
## 333 2010 18 3962
## 334 2010 86 3975
## 335 2010 88 4
## 336 2037 50 3967
## 337 2037 52 19
## 338 2065 42 1
## 339 2065 96 4005
## 340 2080 82 3977
## 341 2080 88 3978
## 342 2087 6 4000
## 343 2087 8 4
## 344 2087 10 1
## 345 2087 13 4
## 346 2087 21 1
## 347 2087 31 4013
## 348 2087 35 5
## 349 2087 38 3
## 350 2087 43 11
## 351 2087 49 3988
## 352 2087 56 2
## 353 2087 59 8
## 354 2087 62 4
## 355 2087 66 7
## 356 2087 68 8
## 357 2087 73 7
## 358 2087 75 3
## 359 2087 81 11
## 360 2087 83 7
## 361 2087 88 1
## 362 2087 94 4012
## 363 2128 43 3985
## 364 2175 53 3933
## 365 2181 8 4002
## 366 2181 14 3976
## 367 2181 20 3980
## 368 2181 49 3988
## 369 2181 52 10
## 370 2181 82 3977
## 371 2192 57 3940
## 372 2196 21 3974
## 373 2196 34 3998
## 374 2196 42 3995
## 375 2196 48 3989
## 376 2196 51 5
## 377 2196 54 11
## 378 2196 56 15
## 379 2196 60 6
## 380 2196 74 4001
## 381 2196 79 8
## 382 2196 81 7
## 383 2196 83 1
## 384 2196 88 7
## 385 2196 96 4005
## 386 2200 14 3976
## 387 2200 37 3986
## 388 2266 23 3981
## 389 2266 57 3940
## 390 2266 82 3977
## 391 2300 88 3978
## 392 2322 46 3985
## 393 2327 22 3979
## 394 2380 51 3954
## 395 2428 23 3981
## 396 2428 39 3985
## 397 2428 58 3943
## 398 2430 55 3939
## 399 2430 59 11
## 400 2430 77 3975
## 401 2430 84 3976
## 402 2437 9 4000
## 403 2437 12 4
## 404 2437 18 3962
## 405 2437 87 3968
## 406 2437 96 4005
## 407 2440 92 4005
## 408 2445 15 3971
## 409 2445 53 3933
## 410 2445 66 3988
## 411 2450 59 3946
## 412 2450 73 4004
## 413 2492 24 3985
## 414 2492 34 3998
## 415 2492 49 3988
## 416 2492 59 3946
## 417 2500 87 3968
## 418 2501 86 3975
## 419 2503 43 3985
## 420 2503 53 3933
## 421 2564 37 3986
## 422 2564 48 3989
## 423 2564 51 5
## 424 2568 65 3991
## 425 2568 85 3972
## 426 2578 30 4014
## 427 2578 35 3997
## 428 2578 51 3954
## 429 2585 79 3972
## 430 2588 54 3933
## 431 2588 68 3987
## 432 2589 50 3967
## 433 2611 80 3974
## 434 2615 15 3971
## 435 2615 36 3996
## 436 2615 68 3987
## 437 2615 75 3989
## 438 2615 89 3984
## 439 2709 10 4000
## 440 2714 57 3940
## 441 2745 48 3989
## 442 2745 51 3
## 443 2745 53 15
## 444 2767 81 3979
## 445 2784 22 3979
## 446 2784 55 3939
## 447 2784 68 3987
## 448 2848 40 3981
## 449 2848 55 3939
## 450 2863 34 3998
## 451 2863 37 9
## 452 2863 49 3988
## 453 2863 59 3946
## 454 2863 63 16
## 455 2863 67 2
## 456 2863 76 3985
## 457 2863 78 9
## 458 2863 84 3976
## 459 2863 88 9
## 460 2863 96 4005
## 461 2907 78 3971
## 462 2954 83 3977
## 463 2972 36 3996
## 464 2975 74 4001
## 465 2982 48 3989
## 466 2982 56 3940
## 467 2982 58 9
## 468 2982 70 3991
## 469 2982 73 3
## 470 2982 84 3976
## 471 2990 52 3947
## 472 2990 66 3988
## 473 2990 85 3972
## 474 2990 88 7
## 475 3044 67 3992
## 476 3048 77 3975
## 477 3056 49 3988
## 478 3056 79 3972
## 479 3282 79 3972
## 480 3285 78 3971
## 481 3315 75 3989
## 482 3315 79 8
## 483 3315 87 3968
## 484 3322 50 3967
## 485 3322 60 3946
## 486 3369 60 3946
## 487 3402 18 3962
## 488 3424 50 3967
## 489 3424 52 19
## 490 3424 74 4001
## 491 3424 77 8
## 492 3424 87 3968
## 493 3424 89 14
## 494 3424 94 4012
## 495 3424 96 3
## 496 3439 20 3980
## 497 3439 23 7
## 498 3439 27 10
## 499 3439 37 3986
## 500 3439 47 3985
## 501 3439 49 7
## 502 3439 52 10
## 503 3439 64 3984
## 504 3439 69 3988
## 505 3439 76 3985
## 506 3439 83 3977
## 507 3498 21 3974
## 508 3498 43 3985
## 509 3498 53 3933
## 510 3498 64 3984
## 511 3667 5 3999
## 512 3667 38 3983
## 513 3866 14 3976
## 514 3866 50 3967
## 515 3866 58 3943
## 516 3924 12 4000
## 517 3924 21 8
## 518 3924 24 6
## 519 3924 30 4014
## 520 3924 42 3995
## 521 3924 44 8
## 522 3924 49 2
## 523 3924 57 3940
## 524 3924 70 3991
## 525 3924 78 3971
## 526 3976 88 3978
## 527 4078 13 3993
## 528 4078 16 7
## 529 4078 22 3979
## 530 4078 48 3989
## 531 4078 53 3933
## 532 4078 59 3946
## 533 4078 70 3991
## 534 4088 89 3984
## 535 4126 36 3996
## 536 4137 59 3946
## 537 4137 75 3989
## 538 4137 93 4012
## 539 4285 83 3977
## 540 4287 76 3985
## 541 4288 55 3939
## 542 4294 48 3989
## 543 4294 53 3933
## 544 4294 55 12
## 545 4294 58 12
## 546 4294 70 3991
## 547 4294 76 3985
## 548 4294 79 8
## 549 4294 81 7
## 550 4294 85 6
## 551 4294 88 7
## 552 4295 60 3946
## 553 4309 5 3999
## 554 4309 25 3981
## 555 4309 32 4012
## 556 4309 36 3
## 557 4309 59 3946
## 558 4310 74 4001
## 559 4313 75 3989
## 560 4315 76 3985
## 561 4358 44 3983
## 562 4358 49 2
## 563 4358 75 3989
## 564 4358 80 11
## 565 4358 83 5
## 566 4358 85 13
## 567 4358 94 4012
## 568 4358 96 3
## 569 4452 87 3968
## 570 4452 93 4012
## 571 4457 9 4000
## 572 4457 21 3974
## 573 4457 34 3998
## 574 4457 55 3939
## 575 4457 57 12
## 576 4457 69 3988
## 577 4457 80 3974
## 578 4457 93 4012
## 579 4458 46 3985
## 580 4458 51 3954
## 581 4458 59 11
## 582 4463 69 3988
## 583 4540 69 3988
## 584 4543 6 4
## 585 4543 14 3976
## 586 4543 19 3982
## 587 4543 23 5
## 588 4543 57 3940
## 589 4543 74 4001
## 590 4543 84 3976
## 591 4544 50 3967
## 592 4544 59 3946
## 593 4559 6 4000
## 594 4559 14 3976
## 595 4559 43 3985
## 596 4559 50 3967
## 597 4559 61 3951
## 598 4559 81 3979
## 599 4731 39 3985
## 600 4731 59 3946
## 601 4731 79 3972
## 602 4734 82 3977
## 603 4760 78 3971
## 604 4760 87 3968
## 605 4800 20 3980
## 606 4800 23 7
## 607 4800 32 4012
## 608 4800 50 3967
## 609 4800 53 18
## 610 4800 55 12
## 611 4800 58 12
## 612 4800 66 3988
## 613 4800 69 7
## 614 4803 8 4002
## 615 4808 42 3995
## 616 4817 64 3984
## 617 4817 72 1
## 618 4817 81 3
## 619 4817 84 7
## 620 4817 87 10
## 621 4817 89 14
## 622 4817 92 4
## 623 4843 50 3967
## 624 4853 21 3974
## 625 4853 26 3981
## 626 4853 42 3995
## 627 4853 48 3989
## 628 4853 53 2
## 629 4853 55 12
## 630 4853 66 7
## 631 4853 68 8
## 632 4853 75 3989
## 633 4853 80 1
## 634 4853 93 5
## 635 4853 98 3
## 636 4857 70 3991
## 637 4857 80 3974
## 638 4889 56 3940
## 639 4889 60 9
## 640 4889 68 3987
## 641 4969 33 4007
## 642 4969 36 2
## 643 4969 53 3933
## 644 4969 87 3968
## 645 5054 83 3977
## 646 5158 9 4000
## 647 5158 43 3985
## 648 5158 54 3933
## 649 5158 59 3946
## 650 5158 78 3971
## 651 5200 51 3954
## 652 5200 54 11
## 653 5200 56 15
## 654 5306 43 3985
## 655 5338 10 1
## 656 5338 12 1
## 657 5338 24 1
## 658 5338 30 4014
## 659 5338 38 3983
## 660 5338 44 3983
## 661 5338 49 3988
## 662 5338 59 3946
## 663 5338 64 14
## 664 5338 74 4001
## 665 5360 87 3968
## 666 5426 22 3979
## 667 5426 52 3947
## 668 5426 56 12
## 669 5431 8 4002
## 670 5431 12 5
## 671 5431 80 3974
## 672 5467 43 3985
## 673 5467 53 3933
## 674 5472 81 3979
## 675 5495 9 4000
## 676 5495 12 4
## 677 5495 19 3982
## 678 5495 23 1
## 679 5495 32 4012
## 680 5495 36 3
## 681 5495 39 10
## 682 5495 43 11
## 683 5495 50 3967
## 684 5495 52 19
## 685 5495 56 2
## 686 5495 59 8
## 687 5495 61 11
## 688 5495 65 4
## 689 5495 67 8
## 690 5495 74 4001
## 691 5495 78 7
## 692 5495 86 3975
## 693 5524 25 3981
## 694 5524 58 3943
## 695 5524 83 3977
## 696 5531 76 3985
## 697 5538 14 3976
## 698 5538 23 3981
## 699 5700 19 3982
## 700 5700 23 5
## 701 5700 46 3985
## 702 5700 53 3933
## 703 5700 57 12
## 704 5702 13 3993
## 705 5714 13 3993
## 706 5714 18 3962
## 707 5714 62 3960
## 708 5763 71 3995
## 709 5844 82 3977
## 710 5899 53 3933
## 711 5947 6 4000
## 712 5947 9 4
## 713 5947 13 1
## 714 5947 17 4
## 715 5947 19 9
## 716 5947 23 3
## 717 5947 28 3991
## 718 5947 40 3981
## 719 5947 43 8
## 720 5947 50 3967
## 721 5947 52 19
## 722 5947 74 4001
## 723 5947 78 7
## 724 5947 96 4005
## 725 6014 9 4000
## 726 6014 14 3976
## 727 6014 20 3980
## 728 6014 26 3981
## 729 6014 36 3996
## 730 6014 43 3985
## 731 6014 51 3954
## 732 6014 55 10
## 733 6014 58 12
## 734 6067 54 3933
## 735 6133 58 3943
## 736 6133 96 4005
## 737 6159 51 3954
## 738 6212 64 3984
## 739 6238 25 3981
## 740 6238 61 3951
## 741 6285 37 3986
## 742 6327 13 3993
## 743 6327 19 3982
## 744 6327 34 3998
## 745 6327 40 3981
## 746 6327 48 3989
## 747 6327 51 3
## 748 6327 54 11
## 749 6327 62 3960
## 750 6327 71 3995
## 751 6327 96 4005
## 752 6354 25 3981
## 753 6354 54 3933
## 754 6354 92 4005
## 755 6444 97 4008
## 756 6447 71 3995
## 757 6462 30 4014
## 758 6462 36 3996
## 759 6462 39 10
## 760 6462 49 3988
## 761 6462 54 3933
## 762 6462 57 13
## 763 6462 63 16
## 764 6462 69 3988
## 765 6462 74 4001
## 766 6462 85 3972
## 767 6462 88 7
## 768 6506 75 3989
## 769 6506 81 3979
## 770 6514 77 3975
## 771 6537 5 3999
## 772 6537 8 4
## 773 6537 43 3985
## 774 6537 84 3976
## 775 6582 13 3993
## 776 6582 39 3985
## 777 6582 52 3947
## 778 6582 57 5
## 779 6582 68 3987
## 780 6582 70 6
## 781 6583 51 3954
## 782 6588 15 3971
## 783 6588 18 7
## 784 6588 21 8
## 785 6588 23 6
## 786 6588 25 3
## 787 6588 39 3985
## 788 6588 44 3983
## 789 6588 48 2
## 790 6588 52 1
## 791 6588 62 3960
## 792 6588 85 3972
## 793 6640 13 3993
## 794 6640 26 3981
## 795 6640 44 3983
## 796 6640 51 3954
## 797 6640 54 11
## 798 6640 56 15
## 799 6640 77 3975
## 800 6641 6 4000
## 801 6641 8 4
## 802 6641 12 5
## 803 6641 16 1
## 804 6641 20 1
## 805 6641 23 7
## 806 6641 25 9
## 807 6641 40 3981
## 808 6641 42 9
## 809 6641 47 3985
## 810 6641 53 3933
## 811 6641 59 2
## 812 6666 73 4004
## 813 6667 14 3976
## 814 6667 52 3947
## 815 6670 17 3964
## 816 6670 26 3981
## 817 6670 60 3946
## 818 6670 77 3975
## 819 6687 62 3960
## 820 6687 66 7
## 821 6687 77 3975
## 822 6687 83 3977
## 823 6743 5 3999
## 824 6802 54 3933
## 825 6819 56 3940
## 826 6840 57 3940
## 827 6937 84 3976
## 828 7024 33 4007
## 829 7024 36 2
## 830 7024 87 3968
## 831 7035 18 3962
## 832 7046 10 4000
## 833 7046 14 2
## 834 7046 33 4007
## 835 7046 39 3985
## 836 7046 50 3967
## 837 7046 54 18
## 838 7046 58 16
## 839 7046 77 3975
## 840 7095 14 3976
## 841 7095 20 3980
## 842 7096 44 3983
## 843 7096 51 3954
## 844 7096 65 3991
## 845 7096 77 3975
## 846 7097 44 3983
## 847 7097 51 3954
## 848 7097 65 3991
## 849 7097 77 3975
## 850 7125 83 3977
## 851 7202 87 3968
## 852 7221 10 4000
## 853 7221 14 2
## 854 7221 33 4007
## 855 7221 39 3985
## 856 7221 50 3967
## 857 7221 54 18
## 858 7221 58 16
## 859 7221 77 3975
## 860 7226 44 3983
## 861 7236 40 3981
## 862 7236 57 3940
## 863 7236 82 3977
## 864 7236 84 9
## 865 7236 89 3984
## 866 7236 93 5
## 867 7249 14 3976
## 868 7314 85 3972
## 869 7345 18 3962
## 870 7345 25 3981
## 871 7345 43 3985
## 872 7345 56 3940
## 873 7367 54 3933
## 874 7401 37 3986
## 875 7404 45 3989
## 876 7449 83 3977
## 877 7568 44 3983
## 878 7568 67 3992
## 879 7583 50 3967
## 880 7601 60 3946
## 881 7601 90 3989
## 882 7798 65 3991
## 883 7798 85 3972
## 884 7860 78 3971
## 885 7876 47 3985
## 886 7897 87 3968
## 887 7902 83 3977
## 888 7946 6 7
## 889 7946 11 4000
## 890 7946 19 8
## 891 7946 21 5
## 892 7946 23 1
## 893 7946 27 10
## 894 7946 29 3
## 895 7946 34 1
## 896 7946 37 9
## 897 7946 47 3985
## 898 7946 49 7
## 899 7946 52 10
## 900 7946 54 11
## 901 7946 72 3996
## 902 7946 85 3972
## 903 7948 90 3989
## 904 8000 14 3976
## 905 8000 29 3996
## 906 8000 42 3995
## 907 8034 11 4000
## 908 8034 14 1
## 909 8034 25 3981
## 910 8034 37 3986
## 911 8034 40 11
## 912 8034 49 3988
## 913 8034 59 3946
## 914 8034 62 4
## 915 8034 71 3995
## 916 8035 35 3997
## 917 8035 60 3946
## 918 8057 33 4007
## 919 8057 36 2
## 920 8057 87 3968
## 921 8060 15 3971
## 922 8060 55 3939
## 923 8257 14 3976
## 924 8288 44 3983
## 925 8288 82 3977
## 926 8288 84 9
## 927 8422 51 3954
## 928 8422 60 3946
## 929 8422 66 7
## 930 8422 68 8
## 931 8422 76 3985
## 932 8422 94 4012
## 933 8432 82 3977
## 934 8502 37 3986
## 935 8512 65 3991
## 936 8513 52 3947
## 937 8513 77 3975
## 938 8519 25 3981
## 939 8519 51 3954
## 940 8519 54 11
## 941 8519 56 15
## 942 8599 17 3964
## 943 8599 19 9
## 944 8599 49 3988
## 945 8599 55 3939
## 946 8637 54 3933
## 947 8694 13 3993
## 948 8694 55 3939
## 949 335 25 0
## 950 410 57 0
## 951 1120 55 0
## 952 1214 19 0
## 953 2065 25 0
## 954 2065 28 0
## 955 2065 30 0
## 956 2065 74 0
## 957 2065 83 0
## 958 2065 86 0
## 959 2065 88 0
## 960 2087 17 0
## 961 2087 25 0
## 962 2087 77 0
## 963 2181 25 0
## 964 2196 85 0
## 965 3924 15 0
## 966 3924 18 0
## 967 4457 25 0
## 968 4458 55 0
## 969 4817 68 0
## 970 4853 60 0
## 971 4853 62 0
## 972 4853 82 0
## 973 5338 7 0
## 974 5338 14 0
## 975 5338 18 0
## 976 5431 84 0
## 977 5947 57 0
## 978 7946 14 0
## 979 8422 79 0
##
## $set.size.summary
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.0 11.5 3971.0 2805.0 3985.0 4014.0
##
## $number.of.treated.units
## [1] 979
##
## $num.units.empty.set
## [1] 31
##
## $lag
## [1] 4#let's examine visually the matched set sizes. The red line shows the number of matched sets without controls or empty matched sets.
plot(msets.maha)
As it looks, most of the sets are of the size between 3500 and 4000.
We can check the covariate balance
get_covariate_balance(PM.results.mah$att,
data = plm_conflict,
use.equal.weights = TRUE,
covariates = c("ELEVATION","POP","ROAD_5KM", "ELECTIONS"),
plot = TRUE)
PE.results <- PanelEstimate(sets = PM.results.mah, data = plm_conflict,
se.method = "bootstrap",
number.iterations = 1000,
confidence.level = .95)
# View the point estimates
PE.results[["estimates"]]## t+0 t+1 t+2 t+3 t+4
## 0.10665310 0.02621308 0.02363874 0.02454290 0.01298473# t+0 t+1 t+2 t+3 t+4
#0.10665310 0.02621308 0.02363874 0.02454290 0.01298473
# View standard errors
PE.results[["standard.error"]]## t+0 t+1 t+2 t+3 t+4
## 0.01843399 0.01463413 0.01482241 0.01601249 0.01416485# t+0 t+1 t+2 t+3 t+4
#0.01898633 0.01408327 0.01502594 0.01623655 0.01473218
#plot the results
plot(PE.results)
As it seems, the effect of the treatment is statistically significant (confidence interval of the treatment [t+0] doesn’t overlap with the control units [t+1,t+n]).
მეჩინგის მთავარი იდეა მარტივია - როცა observational მონაცემები გვაქვს, მანდ მიზეზ-შედეგობრიობაზე საუბარი ძნელია, იმიტომ რომ ვაკვირდებით, ცვლადებს ვერ ვაკონტროლებთ და არ ვერევით გარემოში რომ მანიპულაცია მოვახდინოთ ჩვენთვის სასურველ მოვლენაზე. ამისთვის არსებობს ექსპერიმენტი. მიზეზ-შედეგიბრიობაზე რეალურად ექსპერიმენტის დროს ვლაპარაკობთ. მაგრამ თუ მარტო observational data გვაქვს, შეგვიძლია რაღაცები ვიჩალიჩოთ და მონაცემები დავამსგავსოთ ექსპერიმენტიდან მიღებულ მონაცემებს.
ექსპერიმენტის იდეა რა არის - გყავს ორი იდენტური (თითქმის) ჯგუფი, აქედან ერთს არაფერს უკეთებ (პლაცებო, საკონტროლო ჯგუფი) და მეორეზე რაღაცით მანიპულირებ (წამალს აძლევ, აღქმას უცვლი და ა.შ. ექსპერიმენტული ჯგუფია). თუ გინდა თქვა რომ მანიპულაციამ იმოქმედა შედეგზე, მაშინ უნდა დარწმუნდე რომ სხვა ფაქტორებმა არ იმოქმედა. მაგრამ რახან ერთნაირი ხალხი გვყავს ორივე ჯგუფში, ე.ი. სხვა ფაქტორები გამორიცხულია.
ჩვენი მონაცემები უნდა მივაახლოოთ ამ მდგომარეობასთან. ამისთვის რას ვაკეთებთ - ერთი ცვლადი უნდა ავიღოთ, რომლითაც ვმანიპულირებთ და ვთვლით რომ ეს იწვევს რაღაც მოვლენას (მთავრობის ინტერვენცია - ინსურგენციას). უნდა გავყოთ ორ შესაძლო ვარიანტიად - არსებობს ჩარევა (treatment) ან არ არსებოს (1 და 0 შესაბამისად). ჩვენი მონაცემები მერე ორ ნაწილად შეგვიძლია გავყოთ - ექსპერიმენტული ჯგუფი (treatment ანუ ის დაკვირვებები სადაც 1 გვაქვს) და საკონტროლო ჯგუფი (სადაც არ გვაქვს ინტერვენცია - 0). და ვიღებთ ექსპერიმენტის მსგავს მონაცემებს.
მაგრამ რახან სხვა ფაქტორებს ვერ ვაკონტროლებთ მაგაზეც რაღაც უნდა გავაკეთოთ, რომ მანიპულაციის გარდა სხვა ფაქტორმა არ მოახდინოს გავლენა (გზების მისაწვდომობამ ჩვენ შემთხვევაში მაგალითად). მაგისთვის სპეციალურად ისეთ ქეისებს ვარჩევთ საკონტროლო ჯგუფიდან, რომელიც სხვა ცვლადებზე ძაან გავს ექსპერიმენტულ ჯგუფს. ანუ მაგათ შორის დიდი განსხვავება აღარაა, გარდა იმისა რომ ერთ შემთხვევაში ვმანიპულირებდით, მეორე შემთხვევაში - არა. მაგისთვის ვიყენებთ ლოგისტიკურ რეგრესიას - ვიგებთ თითოეული საკონტროლო ქეისისთვის, რამდენია იმის შანსი რომ 1 იყოს, ანუ იყოს treatment ჯგუფის ნაწილი). ზოგი ქეისი უფრო მიმსგავსებული იქნება ექსპერიმენტულთან, ზოგი ნაკლებად. და ვარჩევთ (პროგრამა არჩევს) საუკეთესო ქეისებს. propensity score ქვია იმ ქულას რითიც ვიგებთ რამდენად მიმსგავსებულია საკონტროლო ქეისებთან.
მერე ადარებ estimate-ს ერთი ჯგუფიდან მეორეს და ნახულობ სტატისტიკური განსხვავება არის თუ არა მაგათ შორის. მაგალითად, ძალადობის საშუალო მაჩვენებელი როგორ განსხვავდება საკონტროლო და ექსპერიმენტულ ჯგუფებში სტატისტიკურად.
მაგრამ ჩვენ შემთხვევაში დროითი მონაცემია (panel data). ამიტომ აქ კიდე ბევრი საჩალიჩოა. ცალკე უნდა გაითვალისწინო რომ ეს ქეისები დროთა განმავლობაში იცვლება და ცვლილებებიც უნდა შეადარო ერთმანეთს. ამისთვის ვიღებთ მაგალითად 4 დროის ცვლილებას (t (იგივე t+0), t+1, t+2, t+3, t+4). ვქმნით სეტებს, სადაც სამი დროის გასვლის დროს ცვლილება არ მომხდარა და მეოთხეზე მოხდა. ანუ, ავისტოში არ ჩაერია მთავრობა, არც სექტემბერში, არც ოქტომბერში, არც ნოემბერში და აგვისტოში - კი. ამ ქეისს ვადარებთ საკონტროლო ქეისებს, სადაც ამ თვეების განმავლობაში მთავრობა საერთოდ არ ჩარეულა. აქამდე გავიგე მეც და მაგის იქით გამიჭირდა ზედაპირული გადაკითხვით და ღრმად აღარ ჩამიხედავს. მაგრამ იდეა სულ ესაა. მერე დეტალებია, ამგალითად რა მეთოდით ამეჩებ ამ ქეისებს. თან ისეა, რომ ერთ ექსპერიმენტულ ქეისს, ბევრი საკონტროლო ქეისი მიესადაგება. და ამათ შორის როგორ უნდა გააკეთო არჩევანი სხვადასხვა მეთოდია. ასევე სხვადასხვა მეთოდია, როგორ აკონტროლებ კოვარიატებს.
სულ ბოლო ფლოთია შენთვის ყველაზე საინტერესო, რაც დასკვნის საშუალებას გაძლევს. მანდ ჩანს რომ ოთხი თვის განმავლობაში რომ არ ჩაერია მთავრობა estimate არის ყველგან 0-ის ტოლი (ნდობის ინტერვალები კვეთს 0-ს, ანუ რეალურად ალბათ 0-ისგან არ განსხვავსება) და როცა მთავრობა ერევა (t+0) ეს estimate 0.1 ან 0.1-ით იმატებს საკონტროლოსგან განსხვავებით. და რახან ნდობის ინტერვალი არ ემთხვევა სხვა დროის estimate-ს ე.ი. განსხვავება სტატისტიკურად მნიშვნელოვანია. ამის ინტერპრეტაცია მაინც გაიკითხე მცოდნე ხალხში. თსუ დამთავრებული ხალხი ვერ გეუბნებით.
წესით, ზემოთ რაც არის კოდი, ეს უნდა დაგჭირდეს. რაღაც ცვლილებების შეტანა შეიძლება მოგიწიოს, მაგალითად lagged ცვლადები ჩასვა კოვარიატებში, ან დროის ინტერვალი აიღო უფრო დიდი ან პატარა. ამაზე ვერაფერს გეტყვი სამწუხაროდ.
თუ ყველაფერი ისეა, როგორც უნდა იყოს, შენი ჰიპოთეზა გამართლდა. ჩვეულებრივ ჰიპოთეზები არ მართლდება ხოლმე, ამიტომ ალბათ რაღაც ვერაა ისე.
მეჩინგისთვის ეს ინსტრუქცია გამოვიყენე: https://cran.r-project.org/web/packages/PanelMatch/vignettes/using_panelmatch.html ამ პაკეტის თეორიული ნაწილი აქედან მოდის: https://web.mit.edu/insong/www/pdf/tscs.pdf