分組運算作法

今天學生問我分組運算的問題,就順便來。舉個最簡單的例子來說,如下表呈現,目標是以山頭(mountain)和坡向(aspect)為一組,來計算thermophilization indicator,計算的方式是把不同的物種依照海拔分布上下限分成不同等級(rank; R),再乘以在該樣區之覆蓋度(cover; C),接下來將所有物種等級覆蓋度乘積加總,再除以該樣區覆蓋度總和,可用下列公式來表示:

\frac{\sum_{i}^{n} C \cdot R}{\sum_{i}^{n} C}

先來簡化問題,若要計算每一座山頭坡向的物種覆蓋度總和,該如何計算。以下用 R 的 data.table 之 group 概念來實作,也就是 DT[, j, by=group_by]:

# 先使用 data.table fread() 來讀取資料
library(data.table)
dummy <- fread('/path/to/the/file.csv', sep = ',')
# 再依照 DT[, j, by=group] 來運算
dummy[, sum(cover), by = .(mountain, aspect)]

就會得到依山頭、方位的覆蓋度總和

mountain aspect V1
1: SEN W11 40.5
2: SEN W13 33.1
3: SUN E11 21.5
4: SUN S13 51.0

因此,我們可以再進一步來計算 thermophilization indicator

Tsi <- dummy[, sum(cover*r)/sum(cover), by = .(mountain,aspect) ]

這個也就是在資料庫結構化查詢語言(Structural Query Language; SQL)的 group by 概念,SQL 寫法如下:

SELECT 
    mountain,
    aspect,
    sum(cover*r)/sum(cover) as tsi
FROM 
    table_name 
GROUP BY 
    mountain,aspect;

範例檔案及 R code:

mountain aspect species cover r
SEN W11 A 30 1
SEN W11 B 10 3
SEN W11 C 0.5 2
SEN W13 A 15 4
SEN W13 D 5 3
SEN W13 E 3 2
SEN W13 F 0.1 1
SEN W13 G 10 4
SUN E11 A 15 2
SUN E11 B 5 3
SUN E11 C 0.5 2
SUN E11 E 1 4
SUN S13 D 35 1
SUN S13 E 10 1
SUN S13 F 5 2
SUN S13 G 1 4
view raw dummy.csv hosted with ❤ by GitHub
library(data.table)
exampleFile <- 'https://gist.githubusercontent.com/mutolisp/4898bb2e25a33487385bf5ede8e553da/raw/f101cc6e126eb2475d2345377bd02bbee3474799/dummy.csv'
dummy <- fread(exampleFile, sep = ',', header = T)
# calculate cover summation by mountain and aspect
dummy[, sum(cover), by = .(mountain, aspect)]
# calculate thermophilization indicator
Tsi <- dummy[, sum(cover*r)/sum(cover), by = .(mountain,aspect) ]
view raw tsi.r hosted with ❤ by GitHub

RMarkdown 轉成 pdf 中文字型問題

Rmarkdown 轉成 pdf 預設會有中文無法顯示的問題(這是因為 knitr 使用 pandoc 轉換成 tex,再用 pdflatex 轉成 pdf)。目前的解法是在 rmarkdown 的 header 加上以下的參數

--
header-includes:
- \usepackage{fontspec} # 使用 fontspec package
- \usepackage{xeCJK}    # 使用 xeCJK package
- \setCJKmainfont{Songti TC} # 指定主要的字型,windows 使用者可用「標楷體」、「新細明體」,或是依照您安裝的字型名稱輸入
output: 
  pdf_document: 
    keep_tex: yes # 保留 tex 檔,萬一出了問題,可以手動檢查並重新編譯
    latex_engine: xelatex # latex 引擎設定為 xelatex
--

之後用 knitr 就會有正常的中文顯示 pdf 了