Tugas 6 | Analisis Regresi dan Prediksi
Komputasi Statistika
| Kontak | : \(\downarrow\) |
| mugemisausan05@gmail.com | |
| https://www.instagram.com/saram.05/ | |
| RPubs | https://rpubs.com/sausanramadhani/ |
Berdasarkan intruksi tugas 6 pada google classroom, diminta untuk mengerjakan tutorial pada dua link yang diberikan (link tertera pada referensi / bagian akhir) di rpubs masing-masing. Berikut ini yang saya pelajari dari kedua tutorial tersebut :
Simple Linier Regression
Dataset berisi observasi mengenai income (antara $15K-$75K) dan happiness (antara 1-10) dan imaginary sample dari 500 orang. Nilai incomedibagi dengan 10.000 untuk mencocokkan data income dengan skala skor happiness (jadi nilai $2 mewakili $20.000, $3 adalah $30.000, dst.)
Getting started in R
library(ggplot2)## Warning: package 'ggplot2' was built under R version 4.1.3library(dplyr)## Warning: package 'dplyr' was built under R version 4.1.3##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(broom)## Warning: package 'broom' was built under R version 4.1.3library(ggpubr)## Warning: package 'ggpubr' was built under R version 4.1.3Step 1: Load the data into R
income.data <- read.csv("income.data.csv")
income.datasummary(income.data)## X income happiness
## Min. : 1.0 Min. :1.506 Min. :0.266
## 1st Qu.:125.2 1st Qu.:3.006 1st Qu.:2.266
## Median :249.5 Median :4.424 Median :3.473
## Mean :249.5 Mean :4.467 Mean :3.393
## 3rd Qu.:373.8 3rd Qu.:5.992 3rd Qu.:4.503
## Max. :498.0 Max. :7.482 Max. :6.863Kedua variabel bersifat kuantitatif. Berdasarkan hasil summary, didapatkan bahwa pada variabel independen (income) :
1. Nilai minimumnya sebesar 1.506
2. Nilai mediannya sebesar 4.424
3. Nilai mean sebesar 4.467
4. Nilai maksimum sebesar 7.482
Sedangkan pada variabel dependen (happiness) :
1. Nilai minimumnya sebesar 0.266
2. Nilai mediannya sebesar 3.473
3. Nilai mean sebesar 3.393
4. Nilai maksimum sebesar 6.863
Step 2: Make sure your data meet the assumptions
Langkah kedua yaitu periksa apakah data tersebut memenuhi empat asumsi utama untuk regresi linier :
1. Independence of Observations (aka no autocorrelation)
Pada dataset income hanya memiliki satu variabel independen dan satu variabel dependen, sehingga tidak perlu menguji hubungan tersembunyi antar variabel.
Catatan : Jika memiliki autokorelasi dalam variabel maka tidak bisa melanjutkan dengan simple linier regression.
2. Normality
Untuk mengecek apakah variabel dependen merupakan distribusi normal, dengan menggunakan fungsi hist() sebagai berikut :
hist(income.data$happiness)
Observasi ditandai dengan berbentuk lonceng yang mana berarti kita bisa lanjutkan ke tahap berikutnya.
3. Linearity
Hubungan antara variabel independen dan dependen harus linear. Kita gunakan scatter plot ntuk melihat apakah distribusi titik data dengan garis lurus.
plot(happiness ~ income, data = income.data)
Berdasarkan hasil di atas, hubungannya linier. Maka bisa dilanjutkan ke tahap berikutnya.
4. Homoscedasticity (aka homogeneity of variance)
Artinya, prediction error tidak berubah secara signifikan selama range prediksi model. Uji asumsi ini bisa dilakukan setelah menyesuaikan model linear.
Step 3: Perform the linear regression analysis
Setelah memeriksa pemenuhan asumsi, maka langkah selanjutnya yaitu melihat apakah ada hubungan linear antara income dan happiness pada survei ini dari 500 orang dengan rentang incomes dari $15k sampai $75k, dimana happiness diukur pada skala 1 hingga 10.
income.happiness.lm <- lm(happiness ~ income, data = income.data)
summary(income.happiness.lm)##
## Call:
## lm(formula = happiness ~ income, data = income.data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.02479 -0.48526 0.04078 0.45898 2.37805
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.20427 0.08884 2.299 0.0219 *
## income 0.71383 0.01854 38.505 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7181 on 496 degrees of freedom
## Multiple R-squared: 0.7493, Adjusted R-squared: 0.7488
## F-statistic: 1483 on 1 and 496 DF, p-value: < 2.2e-16Output di atas berisi tentang rangkuman :
● residu model / residuals (terdapat di awal penjelasan output)
● Koefisien / Coefficients (terdapat estimasi, standar error, statistik t, dan nilai p)
● diagnostik model (tiga baris terakhir yang menunjukkan nilai p untuk mengetahui apakah model sesuai dengan data dengan baik)
Berdasarkan output, dikatakan bahwa ada significant positive relationship antara income dan happiness (p value < 0.0001), dengan 0.713 unit mengalami peningkatan di happiness untuk setiap unit meningkat di income.
Step 4: Check for homoscedasticity
Langkah ini digunakan untuk memastikan bahwa model data sesuai dengan asumsi homoscedascity model linier.
par(mfrow=c(2,2))
plot(income.happiness.lm)
par(mfrow=c(1,1))Berdasarkan Normal Q-Q plot, didapat bahwa residual nyata dari model data membentuk garis satu ke satu yang hampir sempurna dengan residual teoritis dari model yang sempurna. Dengan demikian, dari residual ini bisa dikatakan bahwa model data memenuhi asumsi homoscedasticity.
Step 5: Visualize the results with a graph
Terdapat 4 langkah untuk memvisualisasikan hasil simple linear regression :
1. Plot the data points on a graph
income.graph <- ggplot(income.data, aes(x=income, y=happiness))+
geom_point()
income.graph
2. Add the linear regression line to the plotted data
income.graph <- income.graph + geom_smooth(method="lm", col="black")
income.graph## `geom_smooth()` using formula = 'y ~ x'
3. Add the equation for the regression line.
income.graph <- income.graph +
stat_regline_equation(label.x = 3, label.y = 7)
income.graph## `geom_smooth()` using formula = 'y ~ x'
4. Make the graph ready for publication
income.graph +
theme_bw() +
labs(title = "Reported happiness as a function of income",
x = "Income (x$10,000)",
y = "Happiness score (0 to 10)")## `geom_smooth()` using formula = 'y ~ x'
Step 6: Report your results
Berdasarkan langkah-langkah yang sudah dilakukan, ditemukan bahwa terdapat hubungan yang signifikan antara income dan happiness yaitu p<0.001 dengan peningkatan 0.73 unit happiness untuk setiap peningkatan income sebesar $10,000.
Multiple Linier Regression
Dataset berisi observasi tentang persetase orang biking ke tempat kerja setiap hari, persentase orang smoking, dan persentase heart disease dalam imaginary sample dari 500 kota.
Step 1: Load the data into R
heart.data <- read.csv("heart.data.csv")
heart.dataDalam dataset heart.data, variabel independennya yaitu smoking dan biking, sedangkan variabel dependennya yaitu heart disease.
summary(heart.data)## X biking smoking heart.disease
## Min. : 1.0 Min. : 1.119 Min. : 0.5259 Min. : 0.5519
## 1st Qu.:125.2 1st Qu.:20.205 1st Qu.: 8.2798 1st Qu.: 6.5137
## Median :249.5 Median :35.824 Median :15.8146 Median :10.3853
## Mean :249.5 Mean :37.788 Mean :15.4350 Mean :10.1745
## 3rd Qu.:373.8 3rd Qu.:57.853 3rd Qu.:22.5689 3rd Qu.:13.7240
## Max. :498.0 Max. :74.907 Max. :29.9467 Max. :20.4535Step 2: Make sure your data meet the assumptions
Seperti pada simple linear regression, Multiple linear regression juga perlu memeriksa apakah data sesuai dengan asumsi.
1. Independence of observations (aka no autocorrelation)
cor(heart.data$biking, heart.data$smoking)## [1] 0.01513618Berdasarkan hasil di atas, korelasi antara biking dan smoking hanya 1,5% maka kedua parameter tersebut bisa dimasukkan dalam model.
2. Normality
hist(heart.data$heart.disease)
Distribusinya berbentuk lonceng, maka kita bisa lanjut ke langkah berikutnya.
3. Linearity
plot(heart.disease ~ biking, data=heart.data)
plot(heart.disease ~ smoking, data=heart.data)
Hubungan antara smoking dan heart disease linear walaupun kurang terlihat jelas.
4. Homoscedasticity
Langkah ini dilakukan setelah membuat model.
Step 3: Perform the linear regression analysis
Tingkat biking antara 1 dan 75%, tingkat smoking antara 0.5 dan 30%, dan tingkat heart disease antara 0.5% dan 20.5%. Memeriksa apakah ada hubungan linier antara variabel tersebut :
heart.disease.lm<-lm(heart.disease ~ biking + smoking, data = heart.data)
summary(heart.disease.lm)##
## Call:
## lm(formula = heart.disease ~ biking + smoking, data = heart.data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1789 -0.4463 0.0362 0.4422 1.9331
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 14.984658 0.080137 186.99 <2e-16 ***
## biking -0.200133 0.001366 -146.53 <2e-16 ***
## smoking 0.178334 0.003539 50.39 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.654 on 495 degrees of freedom
## Multiple R-squared: 0.9796, Adjusted R-squared: 0.9795
## F-statistic: 1.19e+04 on 2 and 495 DF, p-value: < 2.2e-16Step 4: Check for homoscedasticity
par(mfrow=c(2,2))
plot(heart.disease.lm)
par(mfrow=c(1,1))Residu tidak menunjukkan bias, maka model sesuai dengan asumsi homoscedasticity.
Step 5: Visualize the results with a graph
Terdapat 7 langkah :
1. Create a new dataframe with the information needed to plot the model
plotting.data<-expand.grid(
biking = seq(min(heart.data$biking), max(heart.data$biking), length.out=30),
smoking=c(min(heart.data$smoking), mean(heart.data$smoking), max(heart.data$smoking)))2. Predict the values of heart disease based on your linear model
plotting.data$predicted.y <- predict.lm(heart.disease.lm, newdata=plotting.data)3. Round the smoking numbers to two decimals
plotting.data$smoking <- round(plotting.data$smoking, digits = 2)4. Change the ‘smoking’ variable into a factor
plotting.data$smoking <- as.factor(plotting.data$smoking)5. Plot the original data
heart.plot <- ggplot(heart.data, aes(x=biking, y=heart.disease)) +
geom_point()
heart.plot
6. Add the regression lines
heart.plot <- heart.plot +
geom_line(data=plotting.data, aes(x=biking, y=predicted.y, color=smoking), size=1.25)## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## i Please use `linewidth` instead.heart.plot
7. Make the graph ready for publication
heart.plot <-
heart.plot +
theme_bw() +
labs(title = "Rates of heart disease (% of population) \n as a function of biking to work and smoking",
x = "Biking to work (% of population)",
y = "Heart disease (% of population)",
color = "Smoking \n (% of population)")
heart.plot
heart.plot + annotate(geom="text", x=30, y=1.75, label=" = 15 + (-0.2*biking) + (0.178*smoking)")
Step 6: Report your results
Terdapat 500 kota dalam survei ini, ditemukan :
● hubungan yang signifikan antara frekuensi biking dan frekuensi heart disease yang mana p<0
● hubungan yang signifikan antara frekuensi smoking dan frekuensi heart disease yang mana p<0.0001
● penurunan sebesar 0.2% pada frekuensi heart disease untuk setiap peningkatan 1% dalam biking
● peningkatan 0.178% pada frekuensi heart disease untuk setiap peningkatan 1% dalam smoking.
Referensi
● https://www.scribbr.com/statistics/linear-regression-in-r/
● https://www.scribbr.com/statistics/multiple-linear-regression/
Sekian pembahasan (pengerjaan) tugas 6 dalam bab Analisis Regresi dan Prediksi pada mata kuliah Statistika Komputasi - Statistika Bisnis Matana University.
Terima Kasih.