自己整理編寫(xiě)的邏輯回歸模板，作為學(xué)習筆記記錄分享。數據集用的是14個(gè)自變量Xi，一個(gè)因變量Y的australian數據集。
1. 測試集和訓練集3、7分組
- australian <- read.csv("australian.csv",as.is = T,sep=",",header=TRUE)
- #讀取行數
- N = length(australian$Y)
- #ind=1的是0.7概率出現的行，ind=2是0.3概率出現的行
- ind=sample(2,N,replace=TRUE,prob=c(0.7,0.3))
- #生成訓練集(這里訓練集和測試集隨機設置為原數據集的70%，30%)
- aus_train <- australian[ind==1,]
- #生成測試集
- aus_test <- australian[ind==2,]
2.生成模型，結果導出
- #生成logis模型，用glm函數
- #用訓練集數據生成logis模型，用glm函數
- #family：每一種響應分布（指數分布族）允許各種關(guān)聯(lián)函數將均值和線(xiàn)性預測器關(guān)聯(lián)起來(lái)。常用的family：binomal(link='logit')--響應變量服從二項分布，連接函數為logit，即logistic回歸
- pre <- glm(Y ~.,family=binomial(link = "logit"),data = aus_train)
- summary(pre)
-
- #測試集的真實(shí)值
- real <- aus_test$Y
- #predict函數可以獲得模型的預測值。這里預測所需的模型對象為pre，預測對象newdata為測試集,預測所需類(lèi)型type選擇response,對響應變量的區間進(jìn)行調整
- predict. <- predict.glm(pre,type='response',newdata=aus_test)
- #按照預測值為1的概率，>0.5的返回1，其余返回0
- predict =ifelse(predict.>0.5,1,0)
- #數據中加入預測值一列
- aus_test$predict = predict
- #導出結果為csv格式
- #write.csv(aus_test,"aus_test.csv")
3.模型檢驗
- ##模型檢驗
- res <- data.frame(real,predict)
- #訓練數據的行數，也就是樣本數量
- n = nrow(aus_train)
- #計算Cox-Snell擬合優(yōu)度
- R2 <- 1-exp((pre$deviance-pre$null.deviance)/n)
- cat("Cox-Snell R2=",R2,"\n")
- #計算Nagelkerke擬合優(yōu)度，我們在最后輸出這個(gè)擬合優(yōu)度值
- R2<-R2/(1-exp((-pre$null.deviance)/n))
- cat("Nagelkerke R2=",R2,"\n")
- ##模型的其他指標
- #residuals(pre) #殘差
- #coefficients(pre) #系數，線(xiàn)性模型的截距項和每個(gè)自變量的斜率，由此得出線(xiàn)性方程表達式?；蛘邔?xiě)為coef(pre)
- #anova(pre) #方差
4.準確率和精度
- true_value=aus_test[,15]
- predict_value=aus_test[,16]
- #計算模型精確度
- error = predict_value-true_value
- accuracy = (nrow(aus_test)-sum(abs(error)))/nrow(aus_test) #精確度--判斷正確的數量占總數的比例
- #計算Precision，Recall和F-measure
- #一般來(lái)說(shuō)，Precision就是檢索出來(lái)的條目（比如：文檔、網(wǎng)頁(yè)等）有多少是準確的，Recall就是所有準確的條目有多少被檢索出來(lái)了
- #和混淆矩陣結合，Precision計算的是所有被檢索到的item（TP+FP）中,"應該被檢索到的item（TP）”占的比例；Recall計算的是所有檢索到的item（TP）占所有"應該被檢索到的item（TP+FN）"的比例。
- precision=sum(true_value & predict_value)/sum(predict_value) #真實(shí)值預測值全為1 / 預測值全為1 --- 提取出的正確信息條數/提取出的信息條數
- recall=sum(predict_value & true_value)/sum(true_value) #真實(shí)值預測值全為1 / 真實(shí)值全為1 --- 提取出的正確信息條數 /樣本中的信息條數
- #P和R指標有時(shí)候會(huì )出現的矛盾的情況，這樣就需要綜合考慮他們，最常見(jiàn)的方法就是F-Measure（又稱(chēng)為F-Score）
- F_measure=2*precision*recall/(precision+recall) #F-Measure是Precision和Recall加權調和平均，是一個(gè)綜合評價(jià)指標
- #輸出以上各結果
- print(accuracy)
- print(precision)
- print(recall)
- print(F_measure)
- #混淆矩陣，顯示結果依次為T(mén)P、FN、FP、TN
- table(true_value,predict_value)
5.ROC曲線(xiàn)的幾個(gè)方法
- #ROC曲線(xiàn)
- # 方法1
- #install.packages("ROCR")
- library(ROCR)
- pred <- prediction(predict.,true_value) #預測值(0.5二分類(lèi)之前的預測值)和真實(shí)值
- performance(pred,'auc')@y.values #AUC值
- perf <- performance(pred,'tpr','fpr')
- plot(perf)
- #方法2
- #install.packages("pROC")
- library(pROC)
- modelroc <- roc(true_value,predict.)
- plot(modelroc, print.auc=TRUE, auc.polygon=TRUE,legacy.axes=TRUE, grid=c(0.1, 0.2),
- grid.col=c("green", "red"), max.auc.polygon=TRUE,
- auc.polygon.col="skyblue", print.thres=TRUE) #畫(huà)出ROC曲線(xiàn)，標出坐標，并標出AUC的值
- #方法3，按ROC定義
- TPR=rep(0,1000)
- FPR=rep(0,1000)
- p=predict.
- for(i in 1:1000)
- {
- p0=i/1000;
- ypred<-1*(p>p0)
- TPR[i]=sum(ypred*true_value)/sum(true_value)
- FPR[i]=sum(ypred*(1-true_value))/sum(1-true_value)
- }
- plot(FPR,TPR,type="l",col=2)
- points(c(0,1),c(0,1),type="l",lty=2)
6.更換測試集和訓練集的選取方式，采用十折交叉驗證
- australian <- read.csv("australian.csv",as.is = T,sep=",",header=TRUE)
- #將australian數據分成隨機十等分
- #install.packages("caret")
- #固定folds函數的分組
- set.seed(7)
- require(caret)
- folds <- createFolds(y=australian$Y,k=10)
-
- #構建for循環(huán)，得10次交叉驗證的測試集精確度、訓練集精確度
-
- max=0
- num=0
-
- for(i in 1:10){
-
- fold_test <- australian[folds[[i]],] #取folds[[i]]作為測試集
- fold_train <- australian[-folds[[i]],] # 剩下的數據作為訓練集
-
- print("***組號***")
-
- fold_pre <- glm(Y ~.,family=binomial(link='logit'),data=fold_train)
- fold_predict <- predict(fold_pre,type='response',newdata=fold_test)
- fold_predict =ifelse(fold_predict>0.5,1,0)
- fold_test$predict = fold_predict
- fold_error = fold_test[,16]-fold_test[,15]
- fold_accuracy = (nrow(fold_test)-sum(abs(fold_error)))/nrow(fold_test)
- print(i)
- print("***測試集精確度***")
- print(fold_accuracy)
- print("***訓練集精確度***")
- fold_predict2 <- predict(fold_pre,type='response',newdata=fold_train)
- fold_predict2 =ifelse(fold_predict2>0.5,1,0)
- fold_train$predict = fold_predict2
- fold_error2 = fold_train[,16]-fold_train[,15]
- fold_accuracy2 = (nrow(fold_train)-sum(abs(fold_error2)))/nrow(fold_train)
- print(fold_accuracy2)
-
-
- if(fold_accuracy>max)
- {
- max=fold_accuracy
- num=i
- }
-
- }
-
- print(max)
- print(num)
-
- ##結果可以看到，精確度accuracy最大的一次為max,取folds[[num]]作為測試集，其余作為訓練集。
7.得到十折交叉驗證的精確度，結果導出
- #十折里測試集最大精確度的結果
- testi <- australian[folds[[num]],]
- traini <- australian[-folds[[num]],] # 剩下的folds作為訓練集
- prei <- glm(Y ~.,family=binomial(link='logit'),data=traini)
- predicti <- predict.glm(prei,type='response',newdata=testi)
- predicti =ifelse(predicti>0.5,1,0)
- testi$predict = predicti
- #write.csv(testi,"ausfold_test.csv")
- errori = testi[,16]-testi[,15]
- accuracyi = (nrow(testi)-sum(abs(errori)))/nrow(testi)
-
- #十折里訓練集的精確度
- predicti2 <- predict.glm(prei,type='response',newdata=traini)
- predicti2 =ifelse(predicti2>0.5,1,0)
- traini$predict = predicti2
- errori2 = traini[,16]-traini[,15]
- accuracyi2 = (nrow(traini)-sum(abs(errori2)))/nrow(traini)
-
- #測試集精確度、取第i組、訓練集精確
- accuracyi;num;accuracyi2
- #write.csv(traini,"ausfold_train.csv")
本站僅提供存儲服務(wù)，所有內容均由用戶(hù)發(fā)布，如發(fā)現有害或侵權內容，請
點(diǎn)擊舉報。
