R ve Kuşlar - eBird Verilerinin R ile Analizi II

R
gis
r-spatial
case_study
R ile Kuş Verilerinin Analizi II
Author

Mehmet Göktuğ Öztürk

Published

December 10, 2023

Merhabalar, blogun ikinci yazısına hoş geldiniz. Takıldığınız ve anlamadığınız yerler olursa lütfen yorum yapmaya çekinmeyiniz. Ayrıca katkılarınızı ve eleştirilerinizi de bekliyorum. Keyifli okumalar.

Blogun ilk yazısını yazarken, ikincisi de birkaç haftaya çıkar, nasıl olsa kodlar hazır diye düşünüyordum ancak tezdir, iştir derken uzadı. Yoğunluk yetmiyormuş gibi bilgisayarın format ihtiyacı ve NTFS dosya sistemine sahip harici diskin linux ile uyumsuzluğu da eklenince yazının taslağı ve kodlarını kaybettim. Sonuç olarak -birbirinden güzel bahanelerimi özetleyecek olursam- gecikmiş, kaybedilmiş ve sıfırdan yazılmış bir yazıyla karşınızdayım.

Bu blog yazısı, bir öncekinin devamı niteliğinde sayılabilir. R ile eBird verilerini işlemeye devam edeceğiz. İlk yazıda birinci soru üzerine odaklanmıştık, bunda ise ikinci soruya odaklanacağız:

  1. Bilindiği gibi Kızılcahamam, Türkiye’deki en önemli kara akbaba - Aegypius monachus popülasyonlarından birisini barındırıyor. Bu türün, Kızılcahamam ilçe sınırları içerisinde nasıl bir dağılımı vardır? Türün dağılımıyla çevresel faktörler arasındaki ilişki kabaca nasıldır?

Temel R ve GIS bilgisine sahip olmanız yazıdaki birçok yeri anlamanız için yeterli olacaktır. Ancak bu alanlarda yeniyseniz de endişelenmeyin. Anlamadığınız her şeyi yorum olarak sorabilirsiniz. Blogun da ana amaçlarından birisi bu zaten. Ayrıca kuşlara ilginiz olmasa bile, bu yazının bazı temel mekânsal analizleri öğrenmek için faydalı olacağını düşünüyorum.

Yazınının Akışı

Yazının genel akışı aşağıdaki gibidir:

  1. Paketlerin yüklenmesi
  2. İlçe verisinin işleri
  3. Kuş verisinin işleri
  4. DEM verisinin işleri
  5. Topografik hesaplamalar
  6. Verinin görselleştirilmesi

1. Paketlerin yüklenmesi

Paket Açıklama
tidyverse Çoğunlukla veri manipülasyonu ve görselleştirme üzerine paketler içeren bir paket koleksiyonu
sf Simple Features: Mekânsal vektör verileri işlemek için
terra Raster verileri işlemek için
tidyterra terra paketi için tidyverse metotları sağlayan bir paket
mapview İnteraktif haritalar için
rgeoboundaries Mülki idare sınırlarını indirmek için

Bu paketler kurulu değilse aşağıdaki kod bloğu ile kurabilirsiniz. Bu kod bloğu, paketi R’a yüklemeye çalışacak, eğer yükleyemezse kuracaktır. Paketler kuruluysa, library() fonksiyonu ile yükleyebiliriz.

# eger paketler yuklu degilse onlari kur
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("sf")) install.packages("sf")
if (!require("terra")) install.packages("terra")
if (!require("tidyterra")) install.packages("tidyterra")
if (!require("mapview")) install.packages("mapview")
if (!require("rgeoboundaries")) remotes::install_github("wmgeolab/rgeoboundaries")

# paketleri R'a yukle
library(tidyverse) # bircok veri isini kolaylastirmak için
library(sf) # r'da mekansal vektor verileri islemek icin
library(terra) # r'da mekansal raster verileri islemek icin
library(tidyterra) # terra ile tidyverse metotlari kullanmak icin
library(mapview) # interaktif haritalar cizmek icin
library(rgeoboundaries) # tr ilce sinirlarina erismek icin

2. İlçe verisinin işleri

İkinci soruya cevap verebilmek için kuş, ilçe sınırları ve dijital yükseklik modeline (DEM) ihtiyacımız var. İlk ilçeden başlayalım. Kuş ve DEM verisini işlemek için ilçeye ihtiyaç duyacağız ne de olsa.

İlk yazıda olduğu gibi bu sefer de geoBoundaries verisini kullanacağız. Şimdi ilçe verilerini indirip, içinden Kızılcahamam’ı çıkartalım.

tr_ilce <- gb_adm2(country = "Turkey", type = "SSCGS") # type = "SSCGS" argumaniyla basitlestirilmiş versiyonunu indiriyoruz
khamam <- tr_ilce[tr_ilce$shapeName == "Kızılcahamam", ]
khamam
Simple feature collection with 1 feature and 6 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 32.27539 ymin: 40.21093 xmax: 32.91882 ymax: 40.75389
Geodetic CRS:  WGS 84
    shapeGroup shapeType    shapeName shapeISO                shapeID
388        TUR      ADM2 Kızılcahamam     <NA> 54988432B2925959455693
    shapeCanonical                       geometry
388      Districts MULTIPOLYGON (((32.51627 40...

Kızılcahamam’ı doğru bir şekilde aldık gibi görünüyor. Verimizi {mapview} paketini kullanarak görselleştirelim, bir kontrol edelim. Bu paket, interaktif haritalar yapmak için oldukça işlevsel.

mapview(khamam)

Gördüğünüz gibi her şey yolunda. :)

3. Kuş verisinin işleri

Kara akbaba. Kaynak: Johann Andreas Naumann’s … Naturgeschichte der Vögel Deutschlands… | Leipzig, G. Fleischer, 1820-1860

Kara akbaba - Aegypius monachus, Avrasya’da geniş dağılıma sahip bir yırtıcı türü. 3 metreye varan kanat açıklığıyla en büyük kuşlardan birisi ve azalan bir popülasyona sahip. Türkiye, türün en önemli yaşam alanlarından birisi. Tür, özellikle Kuzey Anadolu’da, orman - bozkır geçişindeki yaşlı ormanları tercih ediyor. Bu alanlardan birisi de bu yazıda incelediğimiz yer olan Kızılcahamam. Kızılcahamam’ın yaşlı karaçam - Pinus nigra ormanları, türün üremesi için oldukça elverişli.

Türle alakalı kısa bir bilgi verdikten sonra verilere tekrardan dönebiliriz. Kuş verisini nasıl elde ettiğimize ilk yazıda değinmiştik. Tekrardan uzun uzun anlatmaya gerek yok.

Öncelikle kuş verimizi R’a yükleyelim:

ebird <- read_delim("./ebird/ebd_TR_relApr-2023.txt")
print(ebird)
# A tibble: 2,403,720 × 50
   `GLOBAL UNIQUE IDENTIFIER`     `LAST EDITED DATE`  `TAXONOMIC ORDER` CATEGORY
   <chr>                          <dttm>                          <dbl> <chr>   
 1 URN:CornellLabOfOrnithology:E… 2021-04-15 12:59:44             10043 species 
 2 URN:CornellLabOfOrnithology:E… 2021-04-15 12:58:56              5625 species 
 3 URN:CornellLabOfOrnithology:E… 2021-04-15 13:04:32              1361 species 
 4 URN:CornellLabOfOrnithology:E… 2021-04-14 00:03:13             22215 species 
 5 URN:CornellLabOfOrnithology:E… 2021-04-05 11:01:32              5787 species 
 6 URN:CornellLabOfOrnithology:E… 2021-04-21 02:26:13             29239 species 
 7 URN:CornellLabOfOrnithology:E… 2021-04-14 00:03:35             22156 species 
 8 URN:CornellLabOfOrnithology:E… 2021-04-05 11:01:32               291 species 
 9 URN:CornellLabOfOrnithology:E… 2018-09-20 02:46:53              5945 species 
10 URN:CornellLabOfOrnithology:E… 2018-09-20 02:46:53              5945 species 
# ℹ 2,403,710 more rows
# ℹ 46 more variables: `TAXON CONCEPT ID` <chr>, `COMMON NAME` <chr>,
#   `SCIENTIFIC NAME` <chr>, `SUBSPECIES COMMON NAME` <chr>,
#   `SUBSPECIES SCIENTIFIC NAME` <chr>, `EXOTIC CODE` <chr>,
#   `OBSERVATION COUNT` <chr>, `BREEDING CODE` <chr>,
#   `BREEDING CATEGORY` <chr>, `BEHAVIOR CODE` <chr>, `AGE/SEX` <chr>,
#   COUNTRY <chr>, `COUNTRY CODE` <chr>, STATE <chr>, `STATE CODE` <chr>, …

Gördüğünüz gibi 2,403,720 gözlem (satır) ve 50 değişkene (sütun) sahip bir veri tablomuz var. Bu verilerin çoğuna ihtiyaç duymuyoruz. Dolayısıyla işimize yaramayanlardan kurtulmamız işimizi çokça kolaylaştıracaktır. Sadece işimize yarayan sütunları seçip, kara akbaba türüne ait gözlemleri içerecek şekilde filtreleyelim.

karaakbaba <- ebird |>
  select(6, 7, 11, 12, 29, 30) |>  # burada indeks kullanarak sectik, sutun isimleriyle de secebiliriz
  filter(`SCIENTIFIC NAME` == "Aegypius monachus")
karaakbaba
# A tibble: 1,910 × 6
   `COMMON NAME`  `SCIENTIFIC NAME` `OBSERVATION COUNT` `BREEDING CODE` LATITUDE
   <chr>          <chr>             <chr>               <chr>              <dbl>
 1 Cinereous Vul… Aegypius monachus 1                   <NA>                40.5
 2 Cinereous Vul… Aegypius monachus 5                   <NA>                40.5
 3 Cinereous Vul… Aegypius monachus 1                   <NA>                40.0
 4 Cinereous Vul… Aegypius monachus 3                   <NA>                40.5
 5 Cinereous Vul… Aegypius monachus 1                   <NA>                40.5
 6 Cinereous Vul… Aegypius monachus 1                   <NA>                40.5
 7 Cinereous Vul… Aegypius monachus 2                   <NA>                40.5
 8 Cinereous Vul… Aegypius monachus 2                   <NA>                40.5
 9 Cinereous Vul… Aegypius monachus 1                   <NA>                40.5
10 Cinereous Vul… Aegypius monachus 2                   <NA>                40.5
# ℹ 1,900 more rows
# ℹ 1 more variable: LONGITUDE <dbl>

Şimdi değişken sayımız 6’ya, gözlem sayımız ise 1910’a düştü. Veri, çok daha yönetilebilir bir hâlde. Bu etkiyi, değişken isimlerini değiştirerek daha da artırabiliriz. Daha kısa, tamamı küçük harflerden oluşan, latinize ve boşluklar yerine alt çizgi kullandığımız isimler daha kullanışlı olacaktır. Herhangi bir hatayla karşılaşma riskimiz düşecektir. Bunun için ilk olarak karaakbaba data frame’inin değişken (sütun) isimlerine bakalım. names() ya da colnames() fonksiyonunu kullanabiliriz.

names(karaakbaba)
[1] "COMMON NAME"       "SCIENTIFIC NAME"   "OBSERVATION COUNT"
[4] "BREEDING CODE"     "LATITUDE"          "LONGITUDE"

Şimdi de yeni isimlerden oluşan bir karakter vektörü oluşturalım ve bunu karaakbaba data frame’inin değişken isimlerine atayalım.

names(karaakbaba) <- c(
  "eng_name", "sci_name", "obs_count", "breeding_code", "y", "x"
)
names(karaakbaba)
[1] "eng_name"      "sci_name"      "obs_count"     "breeding_code"
[5] "y"             "x"

Gördüğünüz gibi değişken isimlerini istediğimiz gibi değiştirdik.

Artık ilk yüklediğimiz eBird verisini R’dan silebiliriz. Bu veri, 2 milyondan fazla gözlemden oluştuğu için RAM’in şişmesine ve R’ın çökmesine yol açabilir. Ardından da çöp toplayıcıyı (garbage collector) çalıştırmak faydalı olacaktır. R, bu işi genelde otomatik yapıyor ancak büyük veriler silindiği zaman gc()’nin çalıştırılması öneriliyor.

rm(ebird)
gc()
           used  (Mb) gc trigger  (Mb)  max used   (Mb)
Ncells  5214804 278.6   12641914 675.2   9219575  492.4
Vcells 10616247  81.0  127597023 973.5 148118395 1130.1

Şimdi de mekânsallaştırmaya geçebiliriz. Verimiz şu anda bir data frame. Dolayısıyla onunla mekânsal işler yapmak için verimizi mekânsallaştırmamız gerekiyor. Bunun için {sf} paketini kullanacağız. Bu paket, mekânsal vektör veriler işlemek için geliştirmiş, OGC (Open Geospatial Consortium) standardı bir paket. Ben de işlerimin büyük bir kısmında onu tercih ediyorum.

Koordinatların olduğu sütunları ve koordinat sistemini tanımlayarak verimizi sf objesine dönüştürelim.

karaakbaba_sf <- st_as_sf(
  karaakbaba, coords = c("x", "y"), crs = "EPSG:4326"
)
karaakbaba_sf
Simple feature collection with 1910 features and 4 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 26.50652 ymin: 36.02635 xmax: 44.5648 ymax: 42.0362
Geodetic CRS:  WGS 84
# A tibble: 1,910 × 5
   eng_name          sci_name  obs_count breeding_code            geometry
 * <chr>             <chr>     <chr>     <chr>                 <POINT [°]>
 1 Cinereous Vulture Aegypius… 1         <NA>          (32.62995 40.46758)
 2 Cinereous Vulture Aegypius… 5         <NA>          (32.62995 40.46758)
 3 Cinereous Vulture Aegypius… 1         <NA>          (34.61423 40.01914)
 4 Cinereous Vulture Aegypius… 3         <NA>          (32.62995 40.46758)
 5 Cinereous Vulture Aegypius… 1         <NA>           (32.6622 40.48013)
 6 Cinereous Vulture Aegypius… 1         <NA>           (32.62403 40.4557)
 7 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
 8 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
 9 Cinereous Vulture Aegypius… 1         <NA>           (32.62403 40.4557)
10 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
# ℹ 1,900 more rows

Kara akbaba verisi artık mekânsal bir formatta. Geometri tipi, koordinat sistemi gibi öznitelikler barındırıyor.

Veriyi daha iyi tanımak için ilk görselleştirmeleri yapalım. Acaba kara akbabalar Türkiye’de nasıl dağılıyor? Bunun için tekrardan mapview() fonksiyonunu kullanıyoruz.

mapview(karaakbaba_sf)

Bu işlem de tamam. Gördüğünüz gibi ülkemizin kuzeyinde yoğunlaşan bir popülasyon mevcut.

Sıra geldi kuş verimizin Kızılcahamam’a göre alt kümesini almaya.

karaakbaba_subset <- karaakbaba_sf[khamam, ]
karaakbaba_subset
Simple feature collection with 139 features and 4 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 32.29039 ymin: 40.24768 xmax: 32.8661 ymax: 40.64056
Geodetic CRS:  WGS 84
# A tibble: 139 × 5
   eng_name          sci_name  obs_count breeding_code            geometry
   <chr>             <chr>     <chr>     <chr>                 <POINT [°]>
 1 Cinereous Vulture Aegypius… 1         <NA>          (32.62995 40.46758)
 2 Cinereous Vulture Aegypius… 5         <NA>          (32.62995 40.46758)
 3 Cinereous Vulture Aegypius… 3         <NA>          (32.62995 40.46758)
 4 Cinereous Vulture Aegypius… 1         <NA>           (32.6622 40.48013)
 5 Cinereous Vulture Aegypius… 1         <NA>           (32.62403 40.4557)
 6 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
 7 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
 8 Cinereous Vulture Aegypius… 1         <NA>           (32.62403 40.4557)
 9 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
10 Cinereous Vulture Aegypius… 2         <NA>           (32.62403 40.4557)
# ℹ 129 more rows

Kızılcahamam’da, eBird veritabanında kayıtlı 139 tane kara akbaba kaydı varmış. Bunlara bir göz atıp DEM verisiyle ilgilenmeye başlayalım. Kara akbaba dağılımıyla alakalı asıl haritayı en son yapacağız.

plot(
  st_geometry(khamam), 
  main = "Kızılcahamam'da Kara Akbaba Dağılımı", 
  lwd = 2, 
  border = "grey30", 
  axes = TRUE, 
  reset = FALSE
) 
plot(
  st_geometry(karaakbaba_subset), 
  cex = 1, 
  pch = 3, 
  lwd = 1.5, 
  add = TRUE
)

Hemen hızlıca mapview ile de bakalım.

mapview(khamam, lwd = 3, alpha.regions = 0) + mapview(karaakbaba_subset)

4. DEM verisinin işleri

Kara akbabanın Kızılcahamam’daki dağılımını ortaya çıkarttık. Şimdi de bazı çevresel değişkenlerle ilişkisi nasılmış, -kabaca- ona bakalım. Ben, topografik değişkenlerle olan ilişkisini merak ettim ve bu sebeple DEM verisine bakacağız.

DEM (Digital Elevation Model), kısaca, yeryüzünün yükseklik bilgisini sayısal olarak temsil eden modele verilen isim. Her bir grid için, sahip olduğu yükseklik verisini barındırır. Topografik ve hidrolojik birçok analizde sıklıkla kullanılıyor. Rakımın biyolojik sistemlerdeki önemi nedeniyle, ekoloji ve biyocoğrafyada da yaygın bir kullanımı var.

Burada COP-DEM (Copernicus DEM) kullandım. Veriye şu bağlantıdan ulaşabilirsiniz. Belki daha sonra bu verinin nereden, nasıl indirilebileceğiyle alakalı bir başka yazı yazarım.

Raster veri tipi. Kaynak: Geocomputation with R

DEM verisi, raster formatında bir veri. Raster veriyi, bir tür resim ya da bir tür matris olarak düşünebilirsiniz. Bu veri, hücrelerden - gridlerden oluşur ve her bir grid için bir değer içerir. Genellikle, yükseklik, eğim, sıcaklık, yağış gibi sürekli veriler bu formattadır.

Biz bu veri için {terra} paketini kullanacağız. Oldukça hızlı ve işlevsel bir paket.

# dem verisini yukleyelim
dem <- rast("./dem_4326.tif")
dem
class       : SpatRaster 
dimensions  : 7556, 22983, 1  (nrow, ncol, nlyr)
resolution  : 0.0008333333, 0.0008333333  (x, y)
extent      : 25.66542, 44.81792, 35.80792, 42.10458  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : dem_4326.tif 
name        : Copernicus_DSM_COG_30_N35_00_E025_00_DEM 
min value   :                                -86.01817 
max value   :                               5107.41357

Veriyi yazdırdığımızda, çözünürlüğünden, dört köşesinin koordinatına ve koordinat sistemine kadar birçok detayı görebiliyoruz.

plot(dem)

Kolayca haritasını da çizebiliyoruz. Gördüğünüz gibi Türkiye topografik haritası.

Şimdi verimizi Kızılcahamam ilçe sınırlarına göre maskeleyelim. khamam, sf formatında olduğu için vect() ile {terra}’nın kendi vektör veri tipi olan SpatVector’e dönüştürüyoruz.

dem_subset <- crop(dem, vect(khamam), mask = TRUE)
dem_subset
class       : SpatRaster 
dimensions  : 651, 772, 1  (nrow, ncol, nlyr)
resolution  : 0.0008333333, 0.0008333333  (x, y)
extent      : 32.27542, 32.91875, 40.21125, 40.75375  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source(s)   : memory
varname     : dem_4326 
name        : Copernicus_DSM_COG_30_N35_00_E025_00_DEM 
min value   :                                 749.5305 
max value   :                                2069.5696 
plot(dem_subset)

5. Topografik hesaplamalar

terrain() ile eğim ve bakıyı hesaplayalım.

slp <- terrain(dem_subset, "slope")
asp <- terrain(dem_subset, "aspect")
par(mfrow = c(1, 2))
plot(slp, main = "Kızılcahamam Eğim Haritası")
plot(asp, main = "Kızılcahamam Bakı Haritası")

Şimdi de her bir kara akbaba gözlem noktasının sahip olduğu çevresel değişkenleri çıkartalım ve ardından karaakbaba verisiyle birleştirelim.

karaakbaba_subset$ID <- 1:nrow(karaakbaba_subset)

ext_dem <- extract(dem, vect(karaakbaba_subset))
ext_slp <- extract(slp, vect(karaakbaba_subset))
ext_asp <- extract(asp, vect(karaakbaba_subset))

ext_dem <- left_join(karaakbaba_subset, ext_dem, by = "ID")
ext_slp <- left_join(ext_dem, ext_slp, by = "ID")
ext_asp <- left_join(ext_slp, ext_asp, by = "ID")

karaakbaba_ext <- ext_asp
names(karaakbaba_ext)[7] <- "altitude"
glimpse(karaakbaba_ext)
Rows: 139
Columns: 9
$ eng_name      <chr> "Cinereous Vulture", "Cinereous Vulture", "Cinereous Vul…
$ sci_name      <chr> "Aegypius monachus", "Aegypius monachus", "Aegypius mona…
$ obs_count     <chr> "1", "5", "3", "1", "1", "2", "2", "1", "2", "2", "2", "…
$ breeding_code <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ geometry      <POINT [°]> POINT (32.62995 40.46758), POINT (32.62995 40.4675…
$ ID            <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 1…
$ altitude      <dbl> 1058.8070, 1058.8070, 1058.8070, 1045.1802, 1136.2458, 1…
$ slope         <dbl> 11.188035, 11.188035, 11.188035, 22.408733, 11.157285, 1…
$ aspect        <dbl> 309.599090, 309.599090, 309.599090, 242.915925, 7.856463…

Gördüğünüz gibi kara akbaba verimizin içinde artık her bir nokta için rakım, eğim ve bakı bilgisi mevcut. Bu aşamadan sonra modellerle çalışabiliriz ancak bu, bu yazının konusu değil. Ayrıca elimizdeki veri bunun için ne kadar yeterli, ayrı bir tartışma konusu. Şimdi çıkarttığımız topografik verilere bir göz atıp asıl haritamızı yapmaya geçebiliriz.

karaakbaba_ext |>
  st_drop_geometry() |>
  pivot_longer(cols = c("altitude", "slope", "aspect")) |> 
  ggplot(aes(x = value, fill = name, colour = "name")) +
  geom_histogram(alpha = .6) +
  scale_fill_viridis_d() +
  scale_color_viridis_d() +
  theme_minimal() +    
  theme(
    legend.position="none",
  ) + 
  facet_wrap(~name, scales = "free_x")

6. Verinin görselleştirilmesi

Uzun bir yazı oldu yine ancak sonuna geldik. Şimdi Kızılcahamam’daki kara akbaba gözlemlerinin dağılımını, DEM ile birlikte görselleştirmeye çalışalım.

ggplot() +              # grafigi baslatiyor
  geom_spatraster(      # dem verisini ekliyoruz
    data = dem_subset,
    na.rm = TRUE
  ) +
  geom_sf(              # kizilcahamam ilce sinirlarini ekliyoruz
    data = khamam,
    fill = NA
  ) +
  geom_sf(              # kara akbaba gozlem verisini ekliyoruz
    data = karaakbaba_subset,
    aes()
  ) +
  scale_fill_wiki_c() + # dem verisi icin renk paleti ekliyoruz
  labs(                 # etiketleri yaziyoruz
    title = "Kara Akbaba Dağılım Haritası",
    subtitle = "Kızılcahamam'daki kara akbaba gözlemlerinin dağılımı",
    fill = "Yükseklik (m)",
    caption = "Veri: eBird | COP-DEM"
  ) +
  theme_minimal() +     # tema seciyoruz
  theme(                # temanin ozelliklerini belirliyoruz
    plot.background = element_rect("white", colour = "white"),
    text = element_text(family = "ubuntu mono")
  )

Güzel oldu sanki. :)

Bilimle ve huzurla kalınız.