在R中進行線性迴歸時使用實驗不確定度

我有一個我想要擬合多項式的實驗數據集。所述數據包括一個獨立的變量，因變量而在後者的測量的不確定性，例如在R中進行線性迴歸時使用實驗不確定度

2000 0.2084272 0.002067834 
2500 0.207078125 0.001037248 
3000 0.2054202 0.001959138 
3500 0.203488075 0.000328942 
4000 0.2013152 0.000646088 
4500 0.198933825 0.001375657 
5000 0.196375 0.000908696 
5500 0.193668575 0.00014721 
6000 0.1908432 0.000526976 
6500 0.187926325 0.001217318 
7000 0.1849442 0.000556495 
7500 0.181921875 0.000401883 
8000 0.1788832 0.001446992 
8500 0.175850825 0.0
9000 0.1728462 0.001676249 
9500 0.169889575 0.001011735 
10000 0.167  0.000326678

（列X，ÿ，+ -y）。

我可以進行使用上述多項式擬合與例如

mydata = read.table("example.txt") 
model <- lm(V2~V1+I(V1^2)+I(V1^3)+I(V1^4), data = mydata)

但這不使用的不確定性值。我如何告知R，數據集的第三列是不確定性，因此應該用於迴歸分析？

來源

2015-12-22 Joseph Wright

請注意，這裏的數據是'模擬'，基於但不使用真實的東西。 –

您如何使用不確定性？作爲另一個自變量？作爲別的東西？請爲自己做個忙，不要養成附加數據的習慣。它混亂了你的環境，並可能導致分組數據和排序問題。 – Heroka

@Heroka在圖形分析程序（Origin，Igor，...）中有相當標準的分析中使用了一列作爲不確定性：我不是統計學家，所以我不知道它是如何被使用的。在「附加」上，我想你的意思是「附加（數據）」：我從R的書中得知（第二版，_e.g._第467頁），所以假設（d）它是標準的。 –

由於與獨立變量不相關的因變量中的測量誤差，所估計的係數是無偏的，但標準誤差太小。這裏是我使用的參考（第1頁& 2）： http://people.stfx.ca/tleo/econ370term2lec4.pdf

我想你只需要調整由lm（）計算出來的標準錯誤。這就是我在下面的代碼中試圖做的。我不是一個統計人員，所以你可能想要張貼來交叉驗證，並要求更好的直覺。

對於下面的例子，我假定「不確定度」列是標準偏差（或標準誤差）。爲了簡單起見，我將模型改爲：y〜x。

# initialize dataset 
df <- data.frame(
     x = c(2000,2500,3000,3500,4000,4500,5000,5500,6000,6500,7000,7500,8000,8500,9000,9500,10000), 
     y = c(0.2084272,0.207078125,0.2054202,0.203488075,0.2013152,0.198933825, 0.196375,0.193668575, 0.1908432, 0.187926325, 0.1849442, 0.181921875, 0.1788832, 0.175850825, 0.1728462,0.169889575, 0.167), 
     y.err = c(0.002067834, 0.001037248, 0.001959138, 0.000328942, 0.000646088, 0.001375657, 0.000908696, 0.00014721, 0.000526976, 0.001217318, 0.000556495, 0.000401883, 0.001446992, 0.0, 0.001676249, 0.001011735, 0.000326678) 
) 

df 

# model regression 
model <- lm(y~x, data = df) 
summary(model) 

# get the variance of the measurement error 
# thanks to: http://schools-wikipedia.org/wp/v/Variance.htm 
# law of total variance 
pooled.var.y.err <- mean((df$y.err)^2) + var((df$y.err)^2) 

# get variance of beta from model 
# thanks to: http://stats.stackexchange.com/questions/44838/how-are-the-standard-errors-of-coefficients-calculated-in-a-regression 
X <- cbind(1, df$x) 
#  (if you add more variables to the model you need to modify the following line) 
var.betaHat <- anova(model)[[3]][2] * solve(t(X) %*% X) 

# add betaHat variance to measurement error variance 
var.betaHat.adj <- var.betaHat + pooled.var.y.err 

# calculate adjusted standard errors 
sqrt(diag(var.betaHat.adj)) 

# compare to un-adjusted standard errors 
sqrt(diag(var.betaHat))

來源

2016-01-01 19:38:44

在R中進行線性迴歸時使用實驗不確定度

回答

相關問題