Alana Cliantha Dissya - 5003251106

QUESTION 1 : TRIMMED MEAN

A researcher collected monthly household expenditures (in thousands of IDR) from 12 respondents, with 8500 suspected as an outlier. Write an R function trimmed_mean(x, alpha) to compute the ordinary mean (alpha = 0) and the 10% trimmed mean (alpha = 0.10) of the data.

Part (a): Write the trimmed_mean function

trimmed_mean <- function(x, alpha = 0.1) {
  # YOUR CODE HERE
  Y <- sort(x)
  n <- length (x)
  k <- floor(n * alpha)
  
  sumz <-  0
  for (i in (k + 1):(n - k)) {
    sumz <- sumz + Y [i]  
  } 
  return (sumz/(n-2*k))
}

Part (b): Compute ordinary mean and 10% trimmed mean

x <- c(850, 920, 1100, 1250, 1300, 1400, 1500, 1550, 1600, 1750, 2100, 8500)

# Ordinary mean (alpha = 0)
trimmed_mean(x,0)
[1] 1985
# 10% Trimmed mean (alpha = 0.10)

trimmed_mean (x, 0.10)
[1] 1447

DATA VISUALIZATION - QUESTION 1

BOXPLOT

boxplot(x, 
        main = "Box Plot of Monthly Household Expenditures", 
        ylab = "Expenditure (thousands IDR)", 
        col = "pink")

Analysis

The box plot shows that most household expenditures are between 850 and 2100, with one extreme outlier at 8500. The ordinary mean is higher than the majority of the data because of this outlier.

Interpretation

Using a trimmed mean is useful here because it shows a more typical spending value without being messed up by the really high outlier.

HISTOGRAM

hist(x,
     breaks = 15,
     col = "pink",
     border = "black",
     main = "Histogram of Monthly Expenditure",
     xlab = "Expenditure",
     freq = TRUE
)

Analysis

The histogram shows most values are concentrated at lower expenditure levels, with one value far on the right. This indicates the presence of an outlier.

Interpretation

This means the outlier can influence the overall average of the data. Using a trimmed mean helps reduce this effect and gives a more representative result.

QUESTION 2: WEIGHTED MULTIVARIATE DESCRIPTIVE STATISTICS

The Surabaya City Government surveyed households in 31 sub-districts, recording average monthly income, expenditure, and household size. Each sub district is weighted by population. Write an R function weighted_multi_desc(X, w) to calculate weighted descriptive statistics (mean, covariance, and standard deviation) and then apply it to the survey data.

Part (a): Write the weighted_multi_desc function

weighted_multi_desc <- function(X, w) {
  # YOUR CODE HERE
  n_w <-  0
  for (i in 1:length(w)) {
    n_w <- n_w + w [i]
  } 
  
  W = diag(w)
  onevect <- matrix (1, nrow(X), 1)
  
  x_bar_w = (1/n_w) * t(X) %*% W %*% onevect
  
  D <- X - onevect %*% t(x_bar_w)
  
  S_w <- (1/(n_w - 1)) * t(D) %*% W %*% D
  
  sd <- sqrt(diag(S_w))
  
  return(list(W, x_bar_w, S_w, sd))
}

Part (b): Apply to the 31 kecamatan data

# Load data
data <- read.csv("C:/Users/Duta/OneDrive/Documents/Compstat/data_quiz.csv", header = TRUE, sep = ",")
X <- as.matrix(data[, c("Income", "Expenditure", "HH_Size")])
w <- data$Weight


# Apply your function
# YOUR CODE HERE


# Show weighted mean vector
# YOUR CODE HERE
weighted_multi_desc(X, w)[[2]]
                [,1]
Income      4.186869
Expenditure 3.350505
HH_Size     3.832997
# Show weighted covariance matrix
# YOUR CODE HERE
weighted_multi_desc(X, w)[[3]]
                Income Expenditure    HH_Size
Income       0.8703337   0.6840776 -0.4882477
Expenditure  0.6840776   0.5384541 -0.3843073
HH_Size     -0.4882477  -0.3843073  0.2926238
# Show weighted standard deviation vector
# YOUR CODE HERE
weighted_multi_desc(X, w)[[4]]
     Income Expenditure     HH_Size 
  0.9329168   0.7337943   0.5409471 

DATA VISUALIZATION - QUESTION 2

SCATTER PLOT

plot(data$Income, data$Expenditure, 
     main = "Income vs Expenditure (Weighted by Population)",
     xlab = "Income (million IDR)",
     ylab = "Expenditure (million IDR)",
     pch = 19, 
     col = "pink",
     cex = 1.5)

abline(lm(Expenditure ~ Income, data = data),
       col = "blue", lwd = 2)

Analysis

The scatter plot shows an upward trend between income and expenditure. The points generally move in the same direction as income increases.

Interpretation

This means that higher income is associated with higher spending. In general, households with more income tend to spend more.

CORRELATION HEATMAP

install.packages("corrplot")
library(corrplot)

cor_matrix <- cor(data[, c("Income", "Expenditure", "HH_Size")])

corrplot(
  cor_matrix,
  method = "color",        # beda dari "color"
  addCoef.col = "pink",    # beda warna angka
)

Analysis

The plot shows a strong positive relationship between income and expenditure, and weaker relationships with household size.

Interpretation

This means higher income is associated with higher spending, while household size is less related.

LS0tDQp0aXRsZTogIioqUXVpeiAxIERhdGEgVmlzdWFsaXphdGlvbioqIg0Kb3V0cHV0OiBodG1sX25vdGVib29rDQotLS0NCiMjIyMjIF9BbGFuYSBDbGlhbnRoYSBEaXNzeWEgLSA1MDAzMjUxMTA2Xw0KDQojIyBfKipRVUVTVElPTiAxIDogVFJJTU1FRCBNRUFOKipfDQpBIHJlc2VhcmNoZXIgY29sbGVjdGVkIG1vbnRobHkgaG91c2Vob2xkIGV4cGVuZGl0dXJlcyAoaW4gdGhvdXNhbmRzIG9mIElEUikgZnJvbSAxMiByZXNwb25kZW50cywgd2l0aCA4NTAwIHN1c3BlY3RlZCBhcyBhbiBvdXRsaWVyLiBXcml0ZSBhbiBSIGZ1bmN0aW9uIHRyaW1tZWRfbWVhbih4LCBhbHBoYSkgdG8gY29tcHV0ZSB0aGUgb3JkaW5hcnkgbWVhbiAoYWxwaGEgPSAwKSBhbmQgdGhlIDEwJSB0cmltbWVkIG1lYW4gKGFscGhhID0gMC4xMCkgb2YgdGhlIGRhdGEuDQoNCiMjIyAqKlBhcnQgKGEpOiBXcml0ZSB0aGUgdHJpbW1lZF9tZWFuIGZ1bmN0aW9uKioNCmBgYHtyfQ0KdHJpbW1lZF9tZWFuIDwtIGZ1bmN0aW9uKHgsIGFscGhhID0gMC4xKSB7DQogICMgWU9VUiBDT0RFIEhFUkUNCiAgWSA8LSBzb3J0KHgpDQogIG4gPC0gbGVuZ3RoICh4KQ0KICBrIDwtIGZsb29yKG4gKiBhbHBoYSkNCiAgDQogIHN1bXogPC0gIDANCiAgZm9yIChpIGluIChrICsgMSk6KG4gLSBrKSkgew0KICAgIHN1bXogPC0gc3VteiArIFkgW2ldICANCiAgfSANCiAgcmV0dXJuIChzdW16LyhuLTIqaykpDQp9DQoNCmBgYA0KIyMjICoqUGFydCAoYik6IENvbXB1dGUgb3JkaW5hcnkgbWVhbiBhbmQgMTAlIHRyaW1tZWQgbWVhbioqDQpgYGB7cn0NCnggPC0gYyg4NTAsIDkyMCwgMTEwMCwgMTI1MCwgMTMwMCwgMTQwMCwgMTUwMCwgMTU1MCwgMTYwMCwgMTc1MCwgMjEwMCwgODUwMCkNCg0KIyBPcmRpbmFyeSBtZWFuIChhbHBoYSA9IDApDQp0cmltbWVkX21lYW4oeCwwKQ0KDQojIDEwJSBUcmltbWVkIG1lYW4gKGFscGhhID0gMC4xMCkNCg0KdHJpbW1lZF9tZWFuICh4LCAwLjEwKQ0KYGBgDQoNCiMjIF8qKkRBVEEgVklTVUFMSVpBVElPTiAtIFFVRVNUSU9OIDEqKl8NCiMjIyAqKkJPWFBMT1QqKiANCmBgYHtyfQ0KYm94cGxvdCh4LCANCiAgICAgICAgbWFpbiA9ICJCb3ggUGxvdCBvZiBNb250aGx5IEhvdXNlaG9sZCBFeHBlbmRpdHVyZXMiLCANCiAgICAgICAgeWxhYiA9ICJFeHBlbmRpdHVyZSAodGhvdXNhbmRzIElEUikiLCANCiAgICAgICAgY29sID0gInBpbmsiKQ0KYGBgDQojIyMjIyAqKkFuYWx5c2lzKioNClRoZSBib3ggcGxvdCBzaG93cyB0aGF0IG1vc3QgaG91c2Vob2xkIGV4cGVuZGl0dXJlcyBhcmUgYmV0d2VlbiA4NTAgYW5kIDIxMDAsIHdpdGggb25lIGV4dHJlbWUgb3V0bGllciBhdCA4NTAwLiBUaGUgb3JkaW5hcnkgbWVhbiBpcyBoaWdoZXIgdGhhbiB0aGUgbWFqb3JpdHkgb2YgdGhlIGRhdGEgYmVjYXVzZSBvZiB0aGlzIG91dGxpZXIuIA0KDQojIyMjIyAqKkludGVycHJldGF0aW9uKioNClVzaW5nIGEgdHJpbW1lZCBtZWFuIGlzIHVzZWZ1bCBoZXJlIGJlY2F1c2UgaXQgc2hvd3MgYSBtb3JlIHR5cGljYWwgc3BlbmRpbmcgdmFsdWUgd2l0aG91dCBiZWluZyBtZXNzZWQgdXAgYnkgdGhlIHJlYWxseSBoaWdoIG91dGxpZXIuDQoNCg0KIyMjICoqSElTVE9HUkFNKiogDQpgYGB7cn0NCmhpc3QoeCwNCiAgICAgYnJlYWtzID0gMTUsDQogICAgIGNvbCA9ICJwaW5rIiwNCiAgICAgYm9yZGVyID0gImJsYWNrIiwNCiAgICAgbWFpbiA9ICJIaXN0b2dyYW0gb2YgTW9udGhseSBFeHBlbmRpdHVyZSIsDQogICAgIHhsYWIgPSAiRXhwZW5kaXR1cmUiLA0KICAgICBmcmVxID0gVFJVRQ0KKQ0KYGBgDQojIyMjICoqQW5hbHlzaXMqKg0KVGhlIGhpc3RvZ3JhbSBzaG93cyBtb3N0IHZhbHVlcyBhcmUgY29uY2VudHJhdGVkIGF0IGxvd2VyIGV4cGVuZGl0dXJlIGxldmVscywgd2l0aCBvbmUgdmFsdWUgZmFyIG9uIHRoZSByaWdodC4gVGhpcyBpbmRpY2F0ZXMgdGhlIHByZXNlbmNlIG9mIGFuIG91dGxpZXIuDQoNCiMjIyMgKipJbnRlcnByZXRhdGlvbioqDQpUaGlzIG1lYW5zIHRoZSBvdXRsaWVyIGNhbiBpbmZsdWVuY2UgdGhlIG92ZXJhbGwgYXZlcmFnZSBvZiB0aGUgZGF0YS4gVXNpbmcgYSB0cmltbWVkIG1lYW4gaGVscHMgcmVkdWNlIHRoaXMgZWZmZWN0IGFuZCBnaXZlcyBhIG1vcmUgcmVwcmVzZW50YXRpdmUgcmVzdWx0Lg0KDQojIyBfKipRVUVTVElPTiAyOiBXRUlHSFRFRCBNVUxUSVZBUklBVEUgREVTQ1JJUFRJVkUgU1RBVElTVElDUyoqXyANClRoZSBTdXJhYmF5YSBDaXR5IEdvdmVybm1lbnQgc3VydmV5ZWQgaG91c2Vob2xkcyBpbiAzMSBzdWItZGlzdHJpY3RzLCByZWNvcmRpbmcgYXZlcmFnZSBtb250aGx5IGluY29tZSwgZXhwZW5kaXR1cmUsIGFuZCBob3VzZWhvbGQgc2l6ZS4gRWFjaCBzdWIgZGlzdHJpY3QgaXMgd2VpZ2h0ZWQgYnkgcG9wdWxhdGlvbi4gV3JpdGUgYW4gUiBmdW5jdGlvbiB3ZWlnaHRlZF9tdWx0aV9kZXNjKFgsIHcpIHRvIGNhbGN1bGF0ZSB3ZWlnaHRlZCBkZXNjcmlwdGl2ZSBzdGF0aXN0aWNzIChtZWFuLCBjb3ZhcmlhbmNlLCBhbmQgc3RhbmRhcmQgZGV2aWF0aW9uKSBhbmQgdGhlbiBhcHBseSBpdCB0byB0aGUgc3VydmV5IGRhdGEuDQoNCiMjIyAqKlBhcnQgKGEpOiBXcml0ZSB0aGUgd2VpZ2h0ZWRfbXVsdGlfZGVzYyBmdW5jdGlvbioqIA0KYGBge3J9DQp3ZWlnaHRlZF9tdWx0aV9kZXNjIDwtIGZ1bmN0aW9uKFgsIHcpIHsNCiAgIyBZT1VSIENPREUgSEVSRQ0KICBuX3cgPC0gIDANCiAgZm9yIChpIGluIDE6bGVuZ3RoKHcpKSB7DQogICAgbl93IDwtIG5fdyArIHcgW2ldDQogIH0gDQogIA0KICBXID0gZGlhZyh3KQ0KICBvbmV2ZWN0IDwtIG1hdHJpeCAoMSwgbnJvdyhYKSwgMSkNCiAgDQogIHhfYmFyX3cgPSAoMS9uX3cpICogdChYKSAlKiUgVyAlKiUgb25ldmVjdA0KICANCiAgRCA8LSBYIC0gb25ldmVjdCAlKiUgdCh4X2Jhcl93KQ0KICANCiAgU193IDwtICgxLyhuX3cgLSAxKSkgKiB0KEQpICUqJSBXICUqJSBEDQogIA0KICBzZCA8LSBzcXJ0KGRpYWcoU193KSkNCiAgDQogIHJldHVybihsaXN0KFcsIHhfYmFyX3csIFNfdywgc2QpKQ0KfQ0KYGBgDQojIyMgKipQYXJ0IChiKTogQXBwbHkgdG8gdGhlIDMxIGtlY2FtYXRhbiBkYXRhKioNCmBgYHtyfQ0KIyBMb2FkIGRhdGENCmRhdGEgPC0gcmVhZC5jc3YoIkM6L1VzZXJzL0R1dGEvT25lRHJpdmUvRG9jdW1lbnRzL0NvbXBzdGF0L2RhdGFfcXVpei5jc3YiLCBoZWFkZXIgPSBUUlVFLCBzZXAgPSAiLCIpDQpYIDwtIGFzLm1hdHJpeChkYXRhWywgYygiSW5jb21lIiwgIkV4cGVuZGl0dXJlIiwgIkhIX1NpemUiKV0pDQp3IDwtIGRhdGEkV2VpZ2h0DQoNCg0KIyBBcHBseSB5b3VyIGZ1bmN0aW9uDQojIFlPVVIgQ09ERSBIRVJFDQoNCg0KIyBTaG93IHdlaWdodGVkIG1lYW4gdmVjdG9yDQojIFlPVVIgQ09ERSBIRVJFDQp3ZWlnaHRlZF9tdWx0aV9kZXNjKFgsIHcpW1syXV0NCg0KDQojIFNob3cgd2VpZ2h0ZWQgY292YXJpYW5jZSBtYXRyaXgNCiMgWU9VUiBDT0RFIEhFUkUNCndlaWdodGVkX211bHRpX2Rlc2MoWCwgdylbWzNdXQ0KDQoNCiMgU2hvdyB3ZWlnaHRlZCBzdGFuZGFyZCBkZXZpYXRpb24gdmVjdG9yDQojIFlPVVIgQ09ERSBIRVJFDQp3ZWlnaHRlZF9tdWx0aV9kZXNjKFgsIHcpW1s0XV0NCmBgYA0KDQojIyBfKipEQVRBIFZJU1VBTElaQVRJT04gLSBRVUVTVElPTiAyKipfDQojIyMgKipTQ0FUVEVSIFBMT1QqKg0KYGBge3J9DQpwbG90KGRhdGEkSW5jb21lLCBkYXRhJEV4cGVuZGl0dXJlLCANCiAgICAgbWFpbiA9ICJJbmNvbWUgdnMgRXhwZW5kaXR1cmUgKFdlaWdodGVkIGJ5IFBvcHVsYXRpb24pIiwNCiAgICAgeGxhYiA9ICJJbmNvbWUgKG1pbGxpb24gSURSKSIsDQogICAgIHlsYWIgPSAiRXhwZW5kaXR1cmUgKG1pbGxpb24gSURSKSIsDQogICAgIHBjaCA9IDE5LCANCiAgICAgY29sID0gInBpbmsiLA0KICAgICBjZXggPSAxLjUpDQoNCmFibGluZShsbShFeHBlbmRpdHVyZSB+IEluY29tZSwgZGF0YSA9IGRhdGEpLA0KICAgICAgIGNvbCA9ICJibHVlIiwgbHdkID0gMikNCmBgYA0KIyMjIyAqKkFuYWx5c2lzKioNClRoZSBzY2F0dGVyIHBsb3Qgc2hvd3MgYW4gdXB3YXJkIHRyZW5kIGJldHdlZW4gaW5jb21lIGFuZCBleHBlbmRpdHVyZS4gVGhlIHBvaW50cyBnZW5lcmFsbHkgbW92ZSBpbiB0aGUgc2FtZSBkaXJlY3Rpb24gYXMgaW5jb21lIGluY3JlYXNlcy4NCg0KIyMjIyAqKkludGVycHJldGF0aW9uKioNClRoaXMgbWVhbnMgdGhhdCBoaWdoZXIgaW5jb21lIGlzIGFzc29jaWF0ZWQgd2l0aCBoaWdoZXIgc3BlbmRpbmcuIEluIGdlbmVyYWwsIGhvdXNlaG9sZHMgd2l0aCBtb3JlIGluY29tZSB0ZW5kIHRvIHNwZW5kIG1vcmUuDQoNCg0KIyMjICoqQ09SUkVMQVRJT04gSEVBVE1BUCoqDQpgYGB7cn0NCmxpYnJhcnkoY29ycnBsb3QpDQoNCmNvcl9tYXRyaXggPC0gY29yKGRhdGFbLCBjKCJJbmNvbWUiLCAiRXhwZW5kaXR1cmUiLCAiSEhfU2l6ZSIpXSkNCg0KY29ycnBsb3QoDQogIGNvcl9tYXRyaXgsDQogIG1ldGhvZCA9ICJjb2xvciIsDQogIGFkZENvZWYuY29sID0gInBpbmsiLCAgICANCikNCmBgYA0KIyMjIyAqKkFuYWx5c2lzKioNClRoZSBwbG90IHNob3dzIGEgc3Ryb25nIHBvc2l0aXZlIHJlbGF0aW9uc2hpcCBiZXR3ZWVuIGluY29tZSBhbmQgZXhwZW5kaXR1cmUsIGFuZCB3ZWFrZXIgcmVsYXRpb25zaGlwcyB3aXRoIGhvdXNlaG9sZCBzaXplLg0KDQojIyMjICoqSW50ZXJwcmV0YXRpb24qKg0KVGhpcyBtZWFucyBoaWdoZXIgaW5jb21lIGlzIGFzc29jaWF0ZWQgd2l0aCBoaWdoZXIgc3BlbmRpbmcsIHdoaWxlIGhvdXNlaG9sZCBzaXplIGlzIGxlc3MgcmVsYXRlZC4NCg==