R Dersleri - Ali Onur Gitmez

Ana Sayfa Ders 1 Ders 2 Ders 3 Ders 4 Ders 5 Ders 6 Ders 7 Ders 8 Ders 9 Ders 10 Ders 11 Ders 12 Ders 13 Ders 14 Ders 15

Computational Methods

Bu son kisimdaysa gecen hafta causal inference’a giris yaptigimiz gibi sosyal bilimlerde sikca kullanilan computational metodlara bir giris yapip orneklerine bakacagiz. Bu yontemler

  • Web Scraping
  • Text Analysis
  • Sentiment Analysis
  • Network Analysis
  • Spatial Analysis
  • Machine Learning

olacak. Bunlarin hepsi ozellikle son donemlerde yaygin bir sekilde sosyal bilimler alaninda kullaniliyor ve calisilan konuya gore de onemleri yuksek. Bundan dolayı uygulamalı olarak temel bir giriş yapmak önemli olacaktır.

DiÄŸer paketler

library(tidyverse)
library(hms)
library(stringr)
library(writexl)
library(readxl)
library(lubridate)

Scraping paketleri

library(rvest)
library(xml2)
library(XML)
library(httr)

Text paketleri

library(tidytext)
library(stopwords)
library(ggwordcloud)
library(seededlda)
library(glmnet)
library(caret)

Web Scraping

Bir siteden veri cekmek icin once robots.txt kismina bakariz. Orada hangi sayfaların çekmeye uygun olup olmadığını görebiliriz. Bir örnek bile bunu deneyelim.

Sabah Gazetesi’nin sayfasından verileri çekeceğimizi varsayalım. Bunun için öncelikle robots.txt kısmına bakarız. Bunun sonrasında verileri çekebileceğimiz sayfayı bulmamız lazım. Bunun için ya sayfa içinde arama yapmak bir çözüm olabilir. Verileri çekmek için farklı yöntemler kullanabiliriz. Bunlardan ilki HTML ile çekmek olacak. Neredeyse websitesi HTML formatı üstünde bir yapıya kurulu olur ve bunu kullanarak site içinden resim, yazı veya tablo çekebilirsiniz.

h1 ve p gibi kısımlara node deriz. Bunlar veriyi çekmemiz için en önemli kısımlar çünkü bütün veriler bunların altında saklanıyor. Sayfa yapısını inspect ederek veya belli uzantılar aracılığıyla bakarak neyin hangi yapıda olduğunu anlayabiiriz. Şimdi bir örnek ile görelim.

url <- "http://books.toscrape.com/index.html"
webpage <- read_html(url)
titles <- webpage %>%
  html_nodes("h3 a") %>%
  html_attr("title")
titles
##  [1] "A Light in the Attic"                                                                          
##  [2] "Tipping the Velvet"                                                                            
##  [3] "Soumission"                                                                                    
##  [4] "Sharp Objects"                                                                                 
##  [5] "Sapiens: A Brief History of Humankind"                                                         
##  [6] "The Requiem Red"                                                                               
##  [7] "The Dirty Little Secrets of Getting Your Dream Job"                                            
##  [8] "The Coming Woman: A Novel Based on the Life of the Infamous Feminist, Victoria Woodhull"       
##  [9] "The Boys in the Boat: Nine Americans and Their Epic Quest for Gold at the 1936 Berlin Olympics"
## [10] "The Black Maria"                                                                               
## [11] "Starving Hearts (Triangular Trade Trilogy, #1)"                                                
## [12] "Shakespeare's Sonnets"                                                                         
## [13] "Set Me Free"                                                                                   
## [14] "Scott Pilgrim's Precious Little Life (Scott Pilgrim #1)"                                       
## [15] "Rip it Up and Start Again"                                                                     
## [16] "Our Band Could Be Your Life: Scenes from the American Indie Underground, 1981-1991"            
## [17] "Olio"                                                                                          
## [18] "Mesaerion: The Best Science Fiction Stories 1800-1849"                                         
## [19] "Libertarianism for Beginners"                                                                  
## [20] "It's Only the Himalayas"
prices <- webpage %>%
  html_nodes(".price_color") %>%
  html_text() %>%
  sub("£", "", .) %>% 
  as.numeric()   
prices
##  [1] 51.77 53.74 50.10 47.82 54.23 22.65 33.34 17.93 22.60 52.15 13.99 20.66
## [13] 17.46 52.29 35.02 57.25 23.88 37.59 51.33 45.17
ratings <- webpage %>%
  html_nodes(".star-rating") %>%
  html_attr("class") %>%
  sub("star-rating ", "", .)
ratings
##  [1] "Three" "One"   "One"   "Four"  "Five"  "One"   "Four"  "Three" "Four" 
## [10] "One"   "Two"   "Four"  "Five"  "Five"  "Five"  "Three" "One"   "One"  
## [19] "Two"   "Two"
books_df <- data.frame(
  Title = titles,
  Price = prices,
  Rating = ratings,
  stringsAsFactors = FALSE
)

Benzer şekilde websitelerindeki tabloları da çekebiliriz. Burada tabloyu sitede gördüğümüz şekliyle çektiğini görüyoruz.

url2 <- "https://en.wikipedia.org/wiki/List_of_largest_companies_in_the_United_States_by_revenue"
webpage2 <- read_html(url2)
tables <- webpage2 %>% html_nodes("table.wikitable")
companies_table <- tables[[1]] %>% html_table()

Xpath ile de veri çekebiliriz. xpath, html ve xml dökümanları ile çalışmak için güçlü bir tool. Bunu HTML lokasyonlarını öğrenmek için de her objenin kendi yerini öğrenmek için de kullanabiliriz.

url4 <- "http://quotes.toscrape.com/"
webpage4 <- read_html(url4)
quotes <- webpage4 %>%
  html_nodes(xpath = '/html/body/div[1]/div[2]/div[1]/div[1]/span[1]') %>%
  html_text() %>%
  gsub("^\\s+|\\s+$", "", .)

Bu şekilde manüel seçimler yapmak yeriine toplu bir şekilde de veri çekebiliriz. Sabah Gazetesi öğreneğine dönecek olursak:

sitemap_urls <- c("https://www.sabah.com.tr/sitemaparchives/post/2023-1.xml",
                  "https://www.sabah.com.tr/sitemaparchives/post/2023-2.xml",
                  "https://www.sabah.com.tr/sitemaparchives/post/2023-3.xml",
                  "https://www.sabah.com.tr/sitemaparchives/post/2022-12.xml")

2023 yılının ilk 3 ayında Sabah gazetesinde yayınlanmış haberleri çekelim.

df_sabah <- map_dfr(sitemap_urls, function(url) {
  sabah_data <- read_xml(url)
  sabah_urls <- xml_text(xml_find_all(sabah_data, "//d1:loc"))
  sabah_dates <- xml_text(xml_find_all(sabah_data, "//d1:lastmod"))
  sabah_links <- tibble(urls = sabah_urls, dates = sabah_dates)
  sabah_links_2 <- sabah_links %>%
    filter(str_detect(dates, "^2023-(01|02|03)")) %>%
    mutate(dates = substr(dates, 1, 10))
  return(sabah_links_2)
})
num_stories <- nrow(df_sabah)
sprintf("In the first 3 months of 2023, Sabah published %d stories.", num_stories)
## [1] "In the first 3 months of 2023, Sabah published 47228 stories."

Text Analysis

Sonra bu haberlerin kategorizasyonunu yapalım

df_sabah <- df_sabah %>%
  mutate(category = str_extract(urls, "(?<=sabah.com.tr\\/)[^\\/]+")) %>%
  mutate(category = case_when(
    category %in% c("dunya") ~ "World",
    category %in% c("egitim") ~ "Education",
    category %in% c("ekonomi", "finans") ~ "Economy",
    category %in% c("gundem") ~ "Politics",
    category %in% c("kultur-sanat", "magazin", "medya") ~ "Entertainment",
    category %in% c("spor") ~ "Sports",
    category %in% c("yazarlar") ~ "Opinion",
    TRUE ~ "Others"
  ))
df_sabah %>%
  count(category) %>%
  mutate(category = reorder(category, -n)) %>%
  ggplot(aes(x = category, y = n, fill = category)) +
  geom_bar(stat = "identity") +
  ggtitle("Number of News by Category in Sabah Newspaper 2023-01/03") +
  xlab("Category") +
  ylab("Count") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +  
  theme(legend.position = "none")

Haber başlıklarında kaç kere Kılıçdaroğlu kaç kere Erdoğan geçmiş diye bakalım

count_df <- df_sabah %>%
  mutate(names = str_extract_all(urls, "erdogan|kilicdaroglu")) %>%
  unnest(names) %>%
  count(names)
sabah_stories <- df_sabah %>%
  filter(str_detect(urls, "kilicdaroglu|erdogan")) %>%
  distinct(urls, .keep_all = TRUE)
invisible(any(duplicated(sabah_stories$urls)))
paste0("Erdogan appears ", count_df$n[1], 
       " times, and Kilicdaroglu appears ", count_df$n[2], " times.")
## [1] "Erdogan appears 1020 times, and Kilicdaroglu appears 314 times."

Haber içerikleri analizini yapalım

sabah_stories <- sabah_stories %>%
  filter(category == "Politics")
sabah_stories_combined <- data.frame(title = rep(NA, nrow(sabah_stories)),
                                     publish_date = rep(NA, nrow(sabah_stories)),
                                     main_text = rep(NA, nrow(sabah_stories)),
                                     erdogan = rep(NA, nrow(sabah_stories)),
                                     url = rep(NA, nrow(sabah_stories)),
                                     stringsAsFactors = FALSE)
sabah_stories_combined <- sabah_stories_combined %>%
  mutate(publish_date = sabah_stories$dates) %>%
  mutate(urls = sabah_stories$urls) %>%
  mutate(erdogan = ifelse(grepl('erdogan', urls), 1, 0)) %>%
  select(title, publish_date, main_text, erdogan, urls)

Hem haber başlığını hem de textini çıkaralım

for (i in seq_along(sabah_stories$urls)) {
  # get HTML content
  html <- read_html(sabah_stories$urls[i])
  
  # Extract title
  title <- html %>% html_nodes("h1.pageTitle") %>% html_text()
  
  # Extract main text
  main_text_p <- html %>%
    html_nodes("p") %>%
    html_text()
  
  if (length(main_text_p) > 0) {  
    main_text <- main_text_p
  } else { 
    main_text_newsBox <- html %>%
      html_nodes(".newsBox.selectionShareable") %>%
      html_text()
    main_text <- main_text_newsBox
  }
  
  
  sabah_stories_combined$title[i] <- title
  sabah_stories_combined$main_text[i] <- paste(main_text, collapse = "\n")
  
  
  Sys.sleep(7)
}

Kaydedelim ve tekrar çağıralım

write_csv(sabah_stories_combined, "sabahstories.csv")
sabah_stories_combined <- read_csv("sabahstories.csv")
sabah_stories_combined$publish_date <- as.Date(sabah_stories_combined$publish_date)

Günlük olarak çıkan haberlere bakalım

counts <- sabah_stories_combined %>%
  mutate(date = as.Date(publish_date, format = "%Y-%m-%d")) %>%
  group_by(date) %>%
  summarise(erdogan_count = sum(grepl("erdogan", urls, ignore.case = TRUE)),
            kilicdaroglu_count = sum(grepl("kilicdaroglu", urls, ignore.case = TRUE))) %>%
  ungroup()

Başlıklardan noktalamaları çıkaralım, küçük harfe çevirelim ve içinde kılıcdaroğlu geçenleri ayrı bir değişkene alalım.

sabah_stories_combined$title2 <- 
  str_to_lower(str_replace_all(sabah_stories_combined$title, 
                               "[[:punct:]]", ""))
sabah_stories_combined <- sabah_stories_combined %>% 
  mutate(kilicdaroglu = ifelse(grepl('kilicdaroglu', urls), 1, 0))

Stopwords ekleyelim. Bunu R’ın içinde Türkçe stopwords olmadığı için ayrı paketle yapıyoruz.

tr_stopwords <- stopwords::stopwords("tr", source = "stopwords-iso")
tr_stopwords <- tr_stopwords[tr_stopwords != "iyi"]

Sentiment Analysis

Erdoğan geçen başlıklardaki kelimeleri sayalım

erdogan_words <- sabah_stories_combined %>%
  filter(erdogan == 1) %>%
  select(title2) %>%
  unnest_tokens(word, title2) %>%
  filter(!word %in% c("son", "dakika", 
                      tr_stopwords, "recep", "tayyip", "daki̇ka")
         & !grepl("erdoÄŸan", word)) %>%
  count(word, sort = TRUE)

Aynısını Kılıçdaroğlu için yapalım

kilicdaroglu_words <- sabah_stories_combined %>%
  filter(kilicdaroglu == 1) %>%
  select(title2) %>%
  unnest_tokens(word, title2) %>%
  filter(!word %in% c("son", "dakika", 
                      tr_stopwords, "kemal", "daki̇ka") &
           !grepl("kılıçdaroğlu", word)) %>%
  count(word, sort = TRUE)

İkisi için de en çok kullanılan 15 kelimeyi çıkaralım.

top_erdogan_words <- erdogan_words %>%
  slice_head(n = 15)
top_kilicdaroglu_words <- kilicdaroglu_words %>%
  slice_head(n = 15)
erdogan_plot <- ggplot(top_erdogan_words, 
                       aes(x = reorder(word, n), y = n)) +
  geom_bar(stat = "identity", fill = "orange") +
  labs(title = "Most Common Words in Erdogan Titles", 
       y = "Word", x = "Count") +
  theme_minimal() +
  theme(axis.text.y = element_text(angle = 0, hjust = 1)) +
  coord_flip()
ggsave("erdogan.png", dpi = 300)

kilicdaroglu_plot <- ggplot(top_kilicdaroglu_words, 
                            aes(x = reorder(word, n), y = n)) +
  geom_bar(stat = "identity", fill = "red") +
  labs(title = "Most Common Words in Kilicdaroglu Titles", 
       y = "Word", x = "Count") +
  theme_minimal() +
  theme(axis.text.y = element_text(angle = 0, hjust = 1)) +
  coord_flip()
ggsave("kilicdaroglu.png", dpi = 300)

Bunu world cloud olarak da yapabiliriz.

kilicdaroglu_wordcloud <- ggplot(top_kilicdaroglu_words, 
                                 aes(label = word, size = n)) +
  geom_text_wordcloud() +
  scale_size_area(max_size = 15) +
  theme_void() +
  labs(title = "Most Common Words in Kilicdaroglu Titles")
kilicdaroglu_wordcloud

erdogan_wordcloud <- ggplot(top_erdogan_words, 
                                 aes(label = word, size = n)) +
  geom_text_wordcloud() +
  scale_size_area(max_size = 30) +
  theme_void() +
  labs(title = "Most Common Words in Erdogan Titles")
erdogan_wordcloud

Text analizi kısmında textin kendisine ve içinde kelime olarak geçtiğine odaklanmıştık ancak içerik olarak nelerden bahsettiklerine hiç bakmadık.

sentiment <- read_xlsx("sentiment.xlsx")
## New names:
## • `` -> `...1`

Önce verimizi alıp analiz için hazırlayalım

dfm <- dfm(corpus(sentiment$fulltext) %>% 
             tokens(remove_punct = TRUE) %>%
             tokens_tolower() %>%
             tokens_remove(tr_stopwords) %>%
             tokens_wordstem(language = "tr"))

Haberleri konulara bölelim ve beraber bulundukları haberlere göre ortlama olarak nasıl ayrıldıklarına bakalım.

set.seed(94)
lda_model <- textmodel_lda(dfm, k = 10)
top_terms <- terms(lda_model, n= 10)
top_terms_df <- data.frame(top_terms)
top_terms_df_1 <- top_terms_df[ , 1:5]
top_terms_df_2 <- top_terms_df[ , 6:10]
gt_table_1 <- gt::gt(top_terms_df_1) %>% 
  gt::tab_header(title = "Top Terms per Topic (1-5)")
gt_table_2 <- gt::gt(top_terms_df_2) %>% 
  gt::tab_header(title = "Top Terms per Topic (6-10)")
gt_table_1
Top Terms per Topic (1-5)
topic1 topic2 topic3 topic4 topic5
kılıçdaroğlu erdoğa em kemal erdoğa
genel millet erdoğa bay yıl
parti seç çocuk millet ülke
kemal genç dedi işçi millet
başka başka ev iş büyük
chp ak hesap baÅŸka son
kılıçdaroğlu'n ülke medya ifade başka
aday söz atık erdoğa ala
hdp an eÅŸ ada yer
i̇yi̇ parti erdoğan' diyecek dünya
 
gt_table_2
Top Terms per Topic (6-10)
topic6 topic7 topic8 topic9 topic10
başka ülke depre başka su
erdoğa türki erdoğa cumhurbaşkan açılış
recep türk vatandaş terör lira
tayyip dünya başka medya yol
cumhurbaşka uluslararas çalışma erdoğan' başka
genel konu bölge i̇sveç' yüz
ziyaret el konut karar proje
erdoğan' barış felaket örgüt tünel
görüşme halk afet eyle i̇stanbul
kabul ilişki yer i̇sveç yatır

Bu sefer 5 konu altınada toplarlayıp farkına bakalım. Ne değişti?

set.seed(94)
lda_model2 <- textmodel_lda(dfm, k = 5)
top_terms2 <- terms(lda_model2, n = 10)
top_terms2_df <- data.frame(top_terms2)
gt_table_3 <- gt::gt(top_terms2_df) %>% 
  gt::tab_header(title = "Top Terms per Topic")
gt_table_3
Top Terms per Topic
topic1 topic2 topic3 topic4 topic5
baÅŸka erdoÄŸa erdoÄŸa baÅŸka millet
parti depre yıl erdoğa erdoğa
genel başka başka cumhurbaşka ülke
kılıçdaroğlu vatandaş ülke tayyip genç
kemal bölge dünya recep söz
chp çalışma su türk siz
aday konut açılış erdoğan' an
kılıçdaroğlu'n felaket erdoğan' görüşme zama
cumhurbaşka devlet proje türki türki
6 afet tarih kabul dedi
 
topics <- topics(lda_model2)
sentiment$Topic <- topics
topic_labels <-  c("Erdogan","Turkey", "Kilicdaroglu", "Earthquake", "Election" )
sentiment$Topic <- factor(sentiment$Topic, labels = topic_labels)
topic_counts <- table(sentiment$Topic)
topic_counts <- as.data.frame(topic_counts)

Top Terms per Topic
topic1 topic2 topic3 topic4 topic5
an hdp konut emine hizmet
terör masa felaket türk i̇stanbul
i̇sveç akşener depremzede görüşme dünya
14 i̇yi̇ kahramanmaraş dünya alan
bay adaylık ev kabul proje
genç lı afet halk su
mayıs partisi hatay uluslararası kadın
siyasi koalisyon il mesaj sahip
masa isim yardım ilet hayat
örgüt görüşme bura atık yatırım

Konu başına haç haber düştüğünü bu şekilde de görebiliriz

Cluster Frequencies
Cluster Frequency
1 7
2 35
3 6
4 64
5 480

Kelimeleri pozitif ne negatif olarak bölelim. Bunu şimdi öğrenmeyeceğimiz lasso regresyonu diye bir yöntemle yapıyoruz ve böylece hangi kelimelerin neyi tahmin ettiğini görebiliyoruz. Bu kelimeler negatif veya pozitif olarak grupladığımız haberlerde en çok geçen kelimeler.

predictors <- convert(dfm, to = 'matrix')
response <- sentiment$erdogan
set.seed(94)
lasso <- cv.glmnet(predictors, response, family = "binomial", alpha = 1)
tmp_coeffs <- coef(lasso, s = "lambda.min")
coeffs_df <- data.frame(name = tmp_coeffs@Dimnames[[1]][tmp_coeffs@i + 1], coefficient = tmp_coeffs@x)
coeffs_df <- coeffs_df[coeffs_df$name != "(Intercept)",]
coeffs_df <- coeffs_df[order(-abs(coeffs_df$coefficient)),]
top_positive_words <- head(coeffs_df[coeffs_df$coefficient > 0,], 10)
top_negative_words <- head(coeffs_df[coeffs_df$coefficient < 0,], 10)

Negatif haberlerde en çok geçen kelimeler

Top Positive Words
Word Coefficient
ereÄŸli' 7.2696803
beledi̇yeni̇n 2.5043161
baÅŸvura 2.4157332
tazmina 2.3183030
toplanmış 1.7814972
daires 1.0560055
yama 0.9433218
tayyip 0.6628557
etkileÅŸim 0.6190165
onardi 0.4788912

Pozitif haberlerde en çok geçen kelimeler

top_negative_words <- gt::gt(top_negative_words) %>%
  gt::tab_header(title = "Top Negative Words") %>%
  gt::cols_label(name = "Word", coefficient = "Coefficient")
top_negative_words
Top Negative Words
Word Coefficient
kılıçdaroğlu’n -4.749593
sevigen' -4.154531
ağırlayarak -3.702199
flört -3.103029
pervasız -2.786683
kılıçdaroğlu’na -2.731575
tüm -2.563847
mensubiyet -2.469966
çürüdü -2.436991
açiklama -2.348500

İsimleri çıkarıp bir daha bakalım

coeffs_df <- coeffs_df[-grep('erdoğa|kılıçdaroğlu|tayyip|recep|kemal', coeffs_df$name), ]
coeffs_df <- coeffs_df[order(-abs(coeffs_df$coefficient)),]
top_words <- head(coeffs_df, 10)

Bunu yaparken dikkat edilmesi gereken noktayı hatırlatmak gerekirse: robots.txt dışına çıkmak etik olarak görülmez ve zaten bir çok sitede buna yönelik engeller vardır. Sürekli veri çekmeye çalışırsanız da sizi IP üstünden geçici veya kalıcı olarak engelleyebilir. Bu da ileride herhangi bir çalışma yapmamızın önüne geçer.

Bazı sitelerse doğrudan scrape edilmeyi istemezler ancak application programming interface (API) aracılığıyla buna izin verirler. Bunlar ücretli veya ücretsiz olabilir ancak düzenli bir veri çekimine olanak tanırlar.

Daha advanced sentiment analiz icin R’in paketleri var ancak bunlar Ingilizce ile calisiyor. Turkce yapmak icinse Python kullanmamiz lazim. Neyse ki R kullanirken Python cagirip isimizi yapitrabiliyoruz.

Once Python islemlerini yapalim

from transformers import AutoModelForSequenceClassification, AutoTokenizer, pipeline
import pandas as pd
import warnings
from openpyxl import Workbook
warnings.filterwarnings('ignore', category=FutureWarning)
warnings.filterwarnings('ignore', category=UserWarning)
model = AutoModelForSequenceClassification.from_pretrained("savasy/bert-base-turkish-sentiment-cased")
tokenizer = AutoTokenizer.from_pretrained("savasy/bert-base-turkish-sentiment-cased")
sa = pipeline("sentiment-analysis", tokenizer=tokenizer, model=model, device=-1)
input_file = "/Users/alionurgitmez/Desktop/R Dersi/Hafta 15/sabah_PS2.xlsx"
output_file = "after_analysis.xlsx"
df = pd.read_excel(input_file)
print(df.head())
##    Unnamed: 0  ...                                lemmatized_fulltext
## 0           0  ...  BaÅŸkan ErdoÄŸan dan Ali Ä°hsan Destici taziye me...
## 1           1  ...  Emine Erdoğan vatandaş böyle seslen Duyarlı he...
## 2           2  ...  Başkan Erdoğan Yeni Azerbaycan Partisi Gençler...
## 3           3  ...  Son dakika Denizli emekçi kadın buluş Başkan E...
## 4           4  ...  O fotoÄŸraf kahraman konuÅŸ Huriye teyze ErdoÄŸan...
## 
## [5 rows x 8 columns]

Burada sentiment paketini calistiyoruz

sentiment_results = df['title'].apply(lambda x: sa(x)[0])
df['sentiment_label'] = [result['label'] for result in sentiment_results]
df['sentiment_score'] = [result['score'] for result in sentiment_results]
df.to_excel(output_file, index=False)

Burada kaydettigimiz veriyi R’a aliyoruz. Boylece R’a geri donmus olduk. Ikisi de ayni anda calisabiliyor.

## [1] "Average positive sentiment score of Kilicdaroglu articles is: 0.85"
## [1] "Average positive sentiment score of Erdogan articles is: 0.91"

Gorsellestirmesini de yapalim

sabah_stories_positive$date <- as.Date(sabah_stories_positive$date)

average_sentiment2 <- sabah_stories_positive %>%
  mutate(week = as.integer(format(date, "%U"))) %>%
  group_by(erdogan, week) %>%
  summarize(average_sentiment = mean(sentiment_score))

ggplot(average_sentiment2, aes(x = week, y = average_sentiment, color = factor(erdogan))) +
  geom_line() +
  labs(x = "Week", y = "Average Sentiment", color = "Leader") +
  scale_color_manual(values = c("red", "orange"), labels = c("Kilicdaroglu", "Erdogan")) +
  scale_x_continuous(breaks = unique(average_sentiment2$week)) +
  theme_minimal()

Spatial Analysis

Coğrafi analiz özellikle IR alanında ve conflict çalışmalarında gelişmekte olan bir analiz yöntemi. Siyasi ve fiziksel harita veya uydu bilgilerini kullanarak çıkarımlar yapmamıza olanak sağlıyor. Text analizi gibi çok geniş bir alan olduğu için kısa bir giriş yapacağım. Öncelikle bu alanda zorunlu olan paketleri yükleyelim.

library(sf)
library(raster)
library(spData)
library(viridis)
library(scales)

Örnek bir veriyi hazırlayıp bakalım. Bu veri North Carolina sudden infant death syndrome (SIDS) oranlarına bakıyor. Aynı zamanda bu veriyi görselleştirelim. st olarak belirttiğimiz için buna göre görselleştirme yaptı.

nc <- st_read(system.file("shapes/sids.shp", package = "spData"))
## Reading layer `sids' from data source 
##   `/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/spData/shapes/sids.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 100 features and 22 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
## CRS:           NA
plot(nc["SID74"])

Her bölge çevresine 10 km bir güvenli alan oluşturalım.

buffered_nc <- st_buffer(nc, dist = 10000)
plot(buffered_nc, max.plot = 22)

Kesişen alanları bulalım. Bu ilk verimiz ve 10 km koyduğumuz veri arasında kesişen alanları gösteriyor.

intersection_nc <- st_intersection(nc, buffered_nc)
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
plot(intersection_nc)
## Warning: plotting the first 10 out of 44 attributes; use max.plot = 44 to plot
## all

Farklı bir görselleştirme alternatifi olarak da bunu kullanabiliriz.

library(tmap)
tm_shape(nc) +
  tm_polygons(col = "SID74", palette = "Blues", title = "SIDS Rate 1974") +
  tm_layout(legend.outside = TRUE)

Haritalar ile çalışmanın farklı bir örneğini de Türkiye seçimleri üstünden görebiliriz. Bunu ilçe bazlı ve il bazlı seçim sonuçları üstünden yapabiliriz.

library(haven)
library(tidyverse)
library(TRmaps)
ilce_election <- read_dta("Elections_wide.dta")
ilce_election <- ilce_election %>% dplyr::select(adm1_tr, adm2_tr, akp_2018P, chp_2018P, mhp_2018P, 
                                          iyi_2018P, hdp_2018P, sp_2018P)

ilce_election$adm2_tr <- gsub(".*MERKEZ*", "MERKEZ", ilce_election$adm2_tr)

convert_to_lower <- function(str) {
  converted_str <- str_to_sentence(str, locale = "tr")  
  return(converted_str)
}
ilce_election$adm1_tr <- convert_to_lower(ilce_election$adm1_tr)
ilce_election$adm2_tr <- convert_to_lower(ilce_election$adm2_tr)

ilce_election$il_ilce <- paste(ilce_election$adm1_tr, ilce_election$adm2_tr, sep = "-")

ilce_election <- ilce_election %>% dplyr::select(-adm1_tr, -adm2_tr)

tr_ilce <- as.data.frame(st_as_sf(tr_ilce))

tr_ilce <- tr_ilce %>% dplyr::select(il_ilce, geometry, tuik_no, Shape_Area, Shape_Leng)

tr_ilce <- tr_ilce[tr_ilce$il_ilce != "Hakkari-Derecik", ]

ilce_election <- ilce_election %>%
  mutate(il_ilce = ifelse(il_ilce == "Ankara-Kazan", "Ankara-Kahramankazan", il_ilce))
ilce_election <- ilce_election %>%
  mutate(il_ilce = ifelse(il_ilce == "Samsun-19 mayıs", "Samsun-19 Mayıs", il_ilce))
ilce_election <- ilce_election %>%
  mutate(il_ilce = ifelse(il_ilce == "Kırıkkale-Bahşili", "Kırıkkale-Bahşılı", il_ilce))

final_ilce_election <- merge(ilce_election, tr_ilce, by.x = "il_ilce", by.y = "il_ilce")
final_ilce_election$cumhur <- final_ilce_election$akp_2018P + final_ilce_election$mhp_2018P
final_ilce_election$millet <- final_ilce_election$chp_2018P + final_ilce_election$iyi_2018P + final_ilce_election$sp_2018P
final_ilce_election$hdp <- final_ilce_election$hdp_2018P 


final_ilce_election$kazanan <- ifelse(final_ilce_election$hdp > final_ilce_election$millet & final_ilce_election$hdp > final_ilce_election$cumhur, "HDP",
                                      ifelse(final_ilce_election$millet > final_ilce_election$hdp & final_ilce_election$millet > final_ilce_election$cumhur, "MILLET",
                                             ifelse(final_ilce_election$cumhur > final_ilce_election$hdp & final_ilce_election$cumhur > final_ilce_election$millet, "CUMHUR", "ESIT")))
final_ilce_election_sf <- st_as_sf(final_ilce_election)
my_colors <- c("HDP" = "#8f1ca5", "MILLET" = "red", "CUMHUR" = "orange", "CHP" = "firebrick1", 
               "MHP" = "darkred", "IYI" = "cyan", "SAADET" = "maroon", "ESIT" = "darkgrey")

ggplot(final_ilce_election_sf) + 
  geom_sf(aes(fill = kazanan)) +
  scale_fill_manual(values = my_colors) +
  labs(fill = "Kazanan") + 
  theme_void() +
  ggtitle("Tum Ittifaklar Beraber")

Spatial data ile bir analiz öğreneğini inceleyelim. Öncelikle kullanılacak verileri alıp dönüşüm yapıyoruz.

url <- "https://www2.census.gov/geo/tiger/TIGER2019/COUNTY/tl_2019_us_county.zip"
download.file(url, destfile = "us_counties.zip")
unzip("us_counties.zip")
counties <- st_read('/Users/alionurgitmez/Desktop/R Dersi/Hafta 15/tl_2019_us_county.shp')
## Reading layer `tl_2019_us_county' from data source 
##   `/Users/alionurgitmez/Desktop/R Dersi/Hafta 15/tl_2019_us_county.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 3233 features and 17 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -179.2311 ymin: -14.60181 xmax: 179.8597 ymax: 71.43979
## Geodetic CRS:  NAD83
counties_transformed <- st_transform(counties, 2163)
covid_url <- "https://raw.githubusercontent.com/nytimes/covid-19-data/master/us-counties.csv"
covid_data <- read.csv(covid_url)
write_csv(covid_data, "covid_data.csv")
covid_data <- read_csv("covid_data.csv")
## Rows: 2502832 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (2): county, state
## dbl  (3): fips, cases, deaths
## date (1): date
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
covid_latest <- covid_data %>%
  group_by(fips) %>%
  filter(date == max(date)) %>%
  ungroup() %>%
  dplyr::select(fips, cases, deaths)
covid_latest$fips <- as.character(covid_latest$fips)
counties_covid <- counties_transformed %>%
  left_join(covid_latest, by = c("GEOID" = "fips"))
counties_covid <- counties_covid %>%
  mutate(cases_per_100k = (cases / as.numeric(ALAND)) * 100000)

Vakaların çok yüksek olduğu yerlere bakalım

hotspots <- counties_covid %>%
  arrange(desc(cases_per_100k)) %>%
  slice_head(n = 10)
print(hotspots[, c("NAME", "STATEFP", "cases", "cases_per_100k")])
## Simple feature collection with 10 features and 4 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 1930000 ymin: -428834.5 xmax: 2352825 ymax: 140953.9
## Projected CRS: NAD27 / US National Atlas Equal Area
##                    NAME STATEFP  cases cases_per_100k
## 1               Suffolk      25 226017      149.81600
## 2                Hudson      34 179135      149.72731
## 3          Philadelphia      42 317808       91.38143
## 4  District of Columbia      11 143943       90.90732
## 5            Alexandria      51  32647       84.40414
## 6             Arlington      51  46198       68.61124
## 7                 Essex      34 222558       68.15279
## 8         Manassas Park      51   3828       58.54358
## 9                Nassau      36 425386       57.70325
## 10                Union      34 151057       56.75346
##                          geometry
## 1  MULTIPOLYGON (((2322932 130...
## 2  MULTIPOLYGON (((2144411 -12...
## 3  MULTIPOLYGON (((2074643 -24...
## 4  MULTIPOLYGON (((1955476 -40...
## 5  MULTIPOLYGON (((1956985 -41...
## 6  MULTIPOLYGON (((1952303 -40...
## 7  MULTIPOLYGON (((2122970 -12...
## 8  MULTIPOLYGON (((1930000 -42...
## 9  MULTIPOLYGON (((2171046 -10...
## 10 MULTIPOLYGON (((2119662 -14...

Yakındaki yerler üstünden vakaların yüksek olduğu yerleri belirleyelim

high_case_counties <- counties_covid %>%
  filter(cases_per_100k > quantile(cases_per_100k, 0.9, na.rm = TRUE))

neighbors <- st_join(counties_covid, high_case_counties, join = st_touches, left = FALSE)

clusters <- neighbors %>%
  group_by(GEOID.x) %>%
  summarise(cluster_size = n())

Şimdiyse elimizdeki verileri görselleştirelim:

library(viridis)
library(scales)
ggplot(data = counties_covid) +
  geom_sf(aes(fill = cases_per_100k), color = NA) +
  scale_fill_viridis(trans = "log", labels = comma, 
                     name = "Cases per 100,000\n(log scale)") +
  theme_minimal() +
  theme(legend.position = "right") +
  labs(title = "COVID-19 Cases per 100,000 Population by US County",
       subtitle = "Data as of latest available date",
       caption = "Source: New York Times COVID-19 Data")