http://netstat.stat.tku.edu.tw
(User from: 54.162.218.214)
機率  |  統計  |  Data Mining  |  資料處理  |  DEMO  |  Data  |  Links  |  Books  |  相關文件  |   ( ENGLISH )

Home

機率
 計算/查表
 隨機亂數
 Density plots
統計
 基本統計
 卡方檢定
 迴歸分析
 ANOVA
 時間數列
 多變量分析
 品質管制
 無母數方法
Data Mining
 決策樹
 Logistic
 判別分析
 集群分析
 類神經網路
 關聯規則分析
 ALL Methods

R 軟體入門: 0.一個簡單的 R 程式

by Steve Chen ( ) (版權所有, 歡迎散佈)

■ 當我們學習一種新的程式語言時, 最好先觀察這個程式語言的 5 個特色:
  • 1. 常數, 一般變數,陣列變數 (array)的寫法
  • 2. 如何作資料的 Input 與 Output
  • 3. 邏輯判斷 (if)
  • 4. 迴圈 (loop. 例如 for, while, until...)
  • 5. 函數或副程式 (function and/or subroutine) 的寫法
以下這個程式包含了 R 語言這 5 個特色的簡單範例:


# -------------------------------------- program begin
# 以 # 開頭的文字是註解 (comment)
#
# 定義一個計算一組資料的樣本標準差的函數 (function)

my.sd <- function(y)
{
  # y 是一個向量 (vector) 變數, 所以不需用迴圈
  n = length(y)
  kkk = sqrt((sum(y^2) - n * mean(y)^2)/(n - 1))
  return(kkk)
} # function 定義結束


# 假設 mydata.txt 這個檔案只包含 1 列 (row), 2 個數字:
#   0  1.0
#
# 從外部檔案讀入 Normal 分配的兩個參數, 
#     並儲存成 parameters 這個向量 (vetor) 變數的 2 個元素

parameters <- scan("mydata.txt")
mu = parameters[1]      # 此時 mu = 0.0
sigma = parameters[2]   # 此時 sigma = 1.0

# 先創造 2 個向量, 各有 100 及 20000 個元素容量
X <- numeric(100)
results <- numeric(20000)
# Loop 1 (迴圈): 產生 (generate) 20000 組樣本
for(i in 1:20000) 
{
  # Loop 2 (迴圈): 每組樣本中, 隨機產生 100 個觀察值
  for(j in 1:100) 
  {
    # 先隨機產生一個 Uniform(0,1) 的觀察值
    u = runif(1)

    # 若 u <= 0.2, 則從 N(0,10) 產生 1 個觀察值,
    #      並存到 X  這個向量變數的第 j 個元素
    # 若 u >  0.2, 則從 N(0,1) 產生 1 個觀察值, 
    #      並存到 X  這個向量變數的第 j 個元素
    if (u <= 0.2)
    {
       X[j] = rnorm(1,mu,10*sigma)
    } else {
       X[j] = rnorm(1,mu,sigma)
    }
  }
    
  # 呼叫剛剛定義的 function, 計算 X 這組樣本的樣本標準差, 
  # 並且儲存成 results 向量變數的第 i 個元素

  results[i] = my.sd(X)
}

# 畫出所有 20000 個樣本標準差 (每組樣本產生 1 個)的直方圖 (histogram)
hist(results)

# 將 20000 個樣本所算出的 20000 個標準差值存入檔案

write(results,file="output1.dat")

# 同上, 但是輸出存成一個直行 (column), 每個值佔 1 個橫列 (rows)
# 所以 output2.txt 共有 20000 列 (rows) 

write(results,file="output2.txt",ncolumns=1)
# -------------------------------------- program end
將以上程式碼存在一個檔案, 例如, test.R . 然後我們可以使用 R 的批次執行模式 (batch mode) 來跑這個程式 (在 batch mode 中, 圖形部分不會顯現出來):
/usr/local/bin/R BATCH --slave --no-save test.R
執行完畢後, 我們會在同樣目錄下找到一個 test.Rout 的執行記錄檔.

■ 註解 :

以上的隨機樣本產生法, 通常是用來模擬 (simulate) 包含 outliers (極端值)
的樣本. 這樣的機率分配稱為 2 個 Normal 分配組成的 mixture (混合分配):
由 N(0,10) 產生極端值 的機率為 p, 由 N(0,1) 產生正常值的機率為 1 - p.
上面程式中, 假設只有 p = 0.2 的機會可能產生出極端值.