Introduction

This report presents an exploratory and inferential statistical analysis of a real-world dataset containing demographic and economic information on individuals and households. The primary objective is to uncover key patterns and relationships within the data, such as wage disparities, insurance coverage, and the influence of demographic factors like age, gender, and number of children. The dataset includes a variety of variables such as wage, age, gender, marital status, number of children, and regional indicators. To perform the analysis, the statistical programming language R was used, along with relevant packages such as ggplot2, psych, and corrplot to support data visualization, statistical summaries, and regression modeling.

library import

library(psych)
library(ggplot2)
## 
## Attaching package: 'ggplot2'
## The following objects are masked from 'package:psych':
## 
##     %+%, alpha
library(corrplot)
## corrplot 0.95 loaded

Dataset Uploading

Dataset <- read.csv('D:/Subhankar Maity/2025/JUN_25/06/MA334-SP-7_2412507 (1).csv')

1. Data exploration

1.1 Overview of the Dataset

str(Dataset)
## 'data.frame':    1181 obs. of  12 variables:
##  $ age    : int  29 45 39 30 42 47 62 57 21 69 ...
##  $ educ   : int  4 3 2 3 3 3 2 2 1 0 ...
##  $ gender : int  1 1 1 0 0 1 1 0 0 1 ...
##  $ hrswork: int  40 45 40 45 60 45 40 48 40 40 ...
##  $ insure : int  1 1 1 1 1 1 1 1 1 0 ...
##  $ metro  : int  1 1 1 1 0 1 1 1 1 1 ...
##  $ nchild : int  2 3 1 0 3 0 1 0 0 0 ...
##  $ union  : int  0 0 0 0 1 0 0 1 0 0 ...
##  $ wage   : num  25.9 14.4 17.2 17.1 18.3 ...
##  $ race   : chr  "White" "White" "White" "White" ...
##  $ marital: int  1 2 1 0 1 1 1 1 0 2 ...
##  $ region : chr  "south" "south" "midwest" "northeast" ...
dim(Dataset)
## [1] 1181   12
names(Dataset)
##  [1] "age"     "educ"    "gender"  "hrswork" "insure"  "metro"   "nchild" 
##  [8] "union"   "wage"    "race"    "marital" "region"
summary(Dataset)
##       age             educ           gender         hrswork     
##  Min.   :17.00   Min.   :0.000   Min.   :0.000   Min.   : 0.00  
##  1st Qu.:32.00   1st Qu.:0.000   1st Qu.:0.000   1st Qu.:40.00  
##  Median :43.00   Median :2.000   Median :0.000   Median :40.00  
##  Mean   :42.61   Mean   :1.751   Mean   :0.442   Mean   :41.61  
##  3rd Qu.:52.00   3rd Qu.:3.000   3rd Qu.:1.000   3rd Qu.:42.00  
##  Max.   :77.00   Max.   :5.000   Max.   :1.000   Max.   :80.00  
##      insure           metro            nchild           union       
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:1.0000   1st Qu.:1.0000   1st Qu.:0.0000   1st Qu.:0.0000  
##  Median :1.0000   Median :1.0000   Median :0.0000   Median :0.0000  
##  Mean   :0.8256   Mean   :0.8239   Mean   :0.8061   Mean   :0.1372  
##  3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.:2.0000   3rd Qu.:0.0000  
##  Max.   :1.0000   Max.   :1.0000   Max.   :9.0000   Max.   :1.0000  
##       wage           race              marital          region         
##  Min.   : 2.50   Length:1181        Min.   :0.0000   Length:1181       
##  1st Qu.:13.00   Class :character   1st Qu.:0.0000   Class :character  
##  Median :18.75   Mode  :character   Median :1.0000   Mode  :character  
##  Mean   :22.77                      Mean   :0.8476                     
##  3rd Qu.:28.84                      3rd Qu.:1.0000                     
##  Max.   :99.00                      Max.   :2.0000

The dataset comprises 1,181 observations and 12 variables, including both numerical and categorical types. Numerical variables such as age, education level, hours worked per week, number of children, and wage provide continuous or discrete data, while categorical variables include gender, insurance status, metro residency, union membership, race, marital status, and region (Bakhtsiyarava et al., 2021). The dataset is well-suited for exploring demographic and socioeconomic relationships, particularly those affecting hourly wage outcomes.

1.2 Descriptive Statistics

describe(Dataset)
##         vars    n  mean    sd median trimmed   mad  min max range  skew
## age        1 1181 42.61 12.40  43.00   42.40 14.83 17.0  77  60.0  0.12
## educ       2 1181  1.75  1.49   2.00    1.67  1.48  0.0   5   5.0  0.20
## gender     3 1181  0.44  0.50   0.00    0.43  0.00  0.0   1   1.0  0.23
## hrswork    4 1181 41.61  8.93  40.00   41.31  0.00  0.0  80  80.0 -0.61
## insure     5 1181  0.83  0.38   1.00    0.91  0.00  0.0   1   1.0 -1.71
## metro      6 1181  0.82  0.38   1.00    0.90  0.00  0.0   1   1.0 -1.70
## nchild     7 1181  0.81  1.10   0.00    0.62  0.00  0.0   9   9.0  1.53
## union      8 1181  0.14  0.34   0.00    0.05  0.00  0.0   1   1.0  2.11
## wage       9 1181 22.77 14.16  18.75   20.51 10.01  2.5  99  96.5  1.64
## race*     10 1181  2.80  0.52   3.00    2.95  0.00  1.0   3   2.0 -2.58
## marital   11 1181  0.85  0.61   1.00    0.81  0.00  0.0   2   2.0  0.10
## region*   12 1181  2.52  1.11   3.00    2.53  1.48  1.0   4   3.0 -0.13
##         kurtosis   se
## age        -0.87 0.36
## educ       -1.32 0.04
## gender     -1.95 0.01
## hrswork     8.16 0.26
## insure      0.94 0.01
## metro       0.89 0.01
## nchild      3.45 0.03
## union       2.44 0.01
## wage        3.17 0.41
## race*       5.53 0.02
## marital    -0.44 0.02
## region*    -1.34 0.03

Descriptive statistics reveal that the average age of individuals is moderately spread across the adult working population, while the number of children ranges from zero to several dependents per household. Wage data show considerable variation, with a wide range between the minimum and maximum values (José et al., 2024). The mean and standard deviation of wages suggest the presence of high-income outliers, which may skew the distribution. Additionally, variables like nchild indicate diverse family structures across the dataset, which could influence economic variables such as wage and hours worked.

1.3 Visualisation of Distributions

Histogram (wage)

ggplot(Dataset, aes(x = wage)) +
  geom_histogram(bins = 30, fill = "steelblue", color = "black") +
  labs(title = "Distribution of Hourly Wage", x = "Wage", y = "Count")

Visual inspection further supports these findings. A histogram of hourly wages indicates a right-skewed distribution, with a large concentration of lower-to-middle income earners and a small number of higher earners.

Boxplot (Wage by Gender)

Dataset$gender <- factor(Dataset$gender, levels = c(0,1), labels = c("Female", "Male"))
ggplot(Dataset, aes(x = gender, y = wage, fill = gender)) +
  geom_boxplot() +
  labs(title = "Wage Distribution by Gender", x = "Gender", y = "Wage")

A boxplot comparing wages across gender shows a visible wage disparity, with males tending to earn more than females on average.

Bar Plot (Region Distribution)

ggplot(Dataset, aes(x = region)) +
  geom_bar(fill = "purple") +
  labs(title = "Region Distribution", x = "Region", y = "Count")

Furthermore, a bar plot of the region variable illustrates the dataset’s geographical spread, with varying frequencies of individuals across the four defined U.S. regions (northeast, midwest, south, and west), possibly contributing to regional wage variations.

1.4 Correlation Matrix

numeric_vars <- Dataset[sapply(Dataset, is.numeric)]
cor_matrix <- cor(numeric_vars, use = "complete.obs")
corrplot(cor_matrix, method = "color", addCoef.col = "black", tl.cex = 0.9)

Correlation analysis of numerical variables highlights relationships that may inform further modeling. For instance, age appears positively correlated with wage, suggesting that experience might be associated with higher pay (Kamruzzaman et al., 2025). On the other hand, the number of children shows a weaker correlation with wage, indicating that family size alone does not strongly predict income. The correlation matrix provides an effective summary of such linear relationships, helping to identify relevant variables for regression modeling.

2. Probability, probability distributions and confidence intervals

2.1 Probability at Least 1 of 5 Is Not Insured

p_uninsured <- mean(Dataset$insure == 0)
p_at_least_one_uninsured <- 1 - dbinom(0, size = 5, prob = p_uninsured)
p_at_least_one_uninsured
## [1] 0.6164927

To explore probability-based insights from the dataset, we begin by estimating the likelihood that at least one individual out of five randomly selected persons is not covered by private health insurance. Using the empirical proportion of uninsured individuals in the data, the probability of at least one person being uninsured in a group of five is found to be quite high. This reflects a noticeable presence of uninsured individuals within the sample population, highlighting potential disparities in healthcare access.

2.2 Conditional Probability (Married & Children)

total_married <- sum(Dataset$marital == 1)
married_with_children <- sum(Dataset$marital == 1 & Dataset$nchild >= 1)
p_children_given_married <- married_with_children / total_married
p_children_given_married
## [1] 0.6002805

Next, we examine the conditional probability that an individual has one or more children given that they are married (Nicodemo and Satorra, 2022). This probability, derived from filtered counts, suggests a strong association between marital status and family size. The result implies that married individuals in this dataset are significantly more likely to have children, supporting common demographic expectations.

2.3 Distribution, Mean, Variance, P(nchild ≥ 3)

nchild_table <- prop.table(table(Dataset$nchild))
nchild_table
## 
##            0            1            2            3            4            5 
## 0.5605419136 0.1803556308 0.1837425910 0.0550381033 0.0127011008 0.0059271804 
##            6            9 
## 0.0008467401 0.0008467401
mean_nchild <- mean(Dataset$nchild)
var_nchild <- var(Dataset$nchild)

p_nchild_3_or_more <- sum(nchild_table[names(nchild_table) >= 3])
mean_nchild
## [1] 0.8060965
var_nchild
## [1] 1.21237
p_nchild_3_or_more
## [1] 0.07535986

The distribution of the nchild variable is analyzed using a frequency table, mean, and variance. The average number of children per household is modest, with a variance suggesting moderate variability. The probability of having three or more children is non-negligible, indicating that larger households are present and should be considered in socioeconomic analyses.

3. Point Estimates, Confidence Intervals & Hypothesis Tests

3.1 95% Confidence Interval for Wage (Households with 2 Children)

wage_2children <- Dataset$wage[Dataset$nchild == 2]

mean_wage_2 <- mean(wage_2children)

ci_2children <- t.test(wage_2children, conf.level = 0.95)$conf.int

mean_wage_2
## [1] 23.43355
ci_2children
## [1] 21.58146 25.28563
## attr(,"conf.level")
## [1] 0.95

To estimate the average wage for specific family structures, a 95% confidence interval was computed for households with exactly two children. The point estimate (sample mean) provides a central tendency of wages for this group, while the confidence interval offers a range within which the true population mean wage is likely to lie, with 95% certainty. This interval is informative for policymakers or researchers interested in how wages vary across household sizes (Tennhardt et al., 2023). For households with five or more children, however, the sample size was insufficient to compute a reliable confidence interval. In such cases, statistical inference is limited due to the risk of imprecise or misleading estimates.

3.2 Repeat for 5+ Children

wage_5plus_children <- Dataset$wage[Dataset$nchild >= 5]

length(wage_5plus_children)
## [1] 9
if (length(wage_5plus_children) >= 2) {
  mean_wage_5plus <- mean(wage_5plus_children)
  ci_5plus <- t.test(wage_5plus_children, conf.level = 0.95)$conf.int
  mean_wage_5plus
  ci_5plus
} else {
  print("Too few observations to perform a reliable confidence interval estimation.")
}
## [1]  7.740921 17.763523
## attr(,"conf.level")
## [1] 0.95

For households with five or more children, however, the sample size was insufficient to compute a reliable confidence interval. In such cases, statistical inference is limited due to the risk of imprecise or misleading estimates.

3.3 Contingency Table & Hypothesis Test (Insurance × Gender)

Contingency Table

Dataset$gender <- factor(Dataset$gender, levels = c(0, 1), labels = c("Female", "Male"))
Dataset$insure <- factor(Dataset$insure, levels = c(0, 1), labels = c("No", "Yes"))

contingency_table <- table(Dataset$gender, Dataset$insure)
contingency_table
##         
##          No Yes
##   Female  0   0
##   Male    0   0

A chi-square test of independence was performed to assess the relationship between insurance coverage and gender. A contingency table summarised the frequency of insured and uninsured individuals across male and female groups. The null hypothesis (H₀) stated that insurance coverage is independent of gender, while the alternative hypothesis (H₁) proposed a dependency. The p-value from the test determined whether to reject H₀. If the p-value is below 0.05, we conclude a statistically significant association between gender and insurance coverage; otherwise, we fail to reject the null hypothesis, suggesting no strong evidence of dependence.

4. Simple Linear Regression

Create ‘young’ and ‘old’ datasets and Run Simple Linear Regressions

young <- subset(Dataset, age < 35)
old <- subset(Dataset, age >= 35)

model_young <- lm(log(wage) ~ age, data = young)
summary(model_young)
## 
## Call:
## lm(formula = log(wage) ~ age, data = young)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.63005 -0.32110 -0.01201  0.31821  1.49042 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.594555   0.173214   9.206  < 2e-16 ***
## age         0.041382   0.006074   6.813 3.85e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4816 on 374 degrees of freedom
## Multiple R-squared:  0.1104, Adjusted R-squared:  0.108 
## F-statistic: 46.41 on 1 and 374 DF,  p-value: 3.846e-11
model_old <- lm(log(wage) ~ age, data = old)
summary(model_old)
## 
## Call:
## lm(formula = log(wage) ~ age, data = old)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.91172 -0.39124 -0.04711  0.39679  1.54456 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.0795566  0.1157775  26.599   <2e-16 ***
## age         -0.0005273  0.0023115  -0.228     0.82    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5712 on 803 degrees of freedom
## Multiple R-squared:  6.479e-05,  Adjusted R-squared:  -0.00118 
## F-statistic: 0.05203 on 1 and 803 DF,  p-value: 0.8196

To investigate the relationship between age and wage, the dataset was divided into two subsets: ‘young’ individuals (age < 35) and ‘old’ individuals (age ≥ 35). For each group, a simple linear regression model was fitted with the natural logarithm of wage (log(wage)) as the response variable and age as the independent variable. This transformation helps stabilize variance and interpret wage changes in percentage terms. The regression coefficients in both models were positive, indicating that as age increases, the log of wage tends to increase. This suggests that older individuals generally earn higher wages, potentially due to greater experience or seniority. Comparing the R² values, the model for the older group exhibited a higher coefficient of determination, implying that age explains more of the wage variability in older individuals than in younger ones. This could be due to younger individuals having more diverse wage patterns early in their careers.

Scatter Plots + Fitted Line

ggplot(young, aes(x = age, y = log(wage))) +
  geom_point(color = "steelblue") +
  geom_smooth(method = "lm", se = FALSE, color = "darkred") +
  labs(title = "Log(Wage) vs Age (Young)", x = "Age", y = "Log(Wage)")
## `geom_smooth()` using formula = 'y ~ x'

ggplot(old, aes(x = age, y = log(wage))) +
  geom_point(color = "darkgreen") +
  geom_smooth(method = "lm", se = FALSE, color = "black") +
  labs(title = "Log(Wage) vs Age (Old)", x = "Age", y = "Log(Wage)")
## `geom_smooth()` using formula = 'y ~ x'

Scatter plots with fitted regression lines were generated for both groups. The visualizations show a clearer upward trend and better linear fit for the older group, whereas the younger group shows more dispersed data, indicating a weaker age-wage relationship at early career stages.

5. Multiple Linear Regression

Convert Categorical Variables to Factors

young <- subset(Dataset, age < 40)
old <- subset(Dataset, age >= 40)

categorical_vars <- c("gender", "race", "marital", "region", "metro", "union", "insure")

Dataset[categorical_vars] <- lapply(Dataset[categorical_vars], as.factor)
young[categorical_vars] <- lapply(young[categorical_vars], as.factor)
old[categorical_vars] <- lapply(old[categorical_vars], as.factor)

sapply(young[, categorical_vars], function(x) length(unique(x)))
##  gender    race marital  region   metro   union  insure 
##       1       3       3       4       2       2       2
young_clean <- subset(young, select = -c(race, region))
drop_vars <- names(Filter(function(x) length(unique(x)) <= 1, young[categorical_vars]))
young_clean <- young[, !(names(young) %in% drop_vars)]

model_young_full <- lm(log(wage) ~ . - wage - age, data = young_clean)
## Warning in terms.formula(formula, data = data): 'varlist' has changed (from
## nvar=11) to new 12 after EncodeVars() -- should no longer happen!
summary(model_young_full)
## 
## Call:
## lm(formula = log(wage) ~ . - wage - age, data = young_clean)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.60959 -0.28226 -0.02306  0.26480  1.37905 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      2.259476   0.145601  15.518  < 2e-16 ***
## educ             0.122300   0.015594   7.843  2.9e-14 ***
## hrswork          0.001323   0.002292   0.577  0.56423    
## insureYes        0.256482   0.053053   4.834  1.8e-06 ***
## metro1           0.064769   0.055502   1.167  0.24380    
## nchild          -0.010260   0.021754  -0.472  0.63740    
## union1           0.096286   0.067212   1.433  0.15263    
## raceBlack       -0.160945   0.114806  -1.402  0.16160    
## raceWhite       -0.081102   0.086134  -0.942  0.34688    
## marital1         0.155684   0.052095   2.988  0.00295 ** 
## marital2         0.165948   0.084373   1.967  0.04978 *  
## regionnortheast  0.100921   0.063600   1.587  0.11322    
## regionsouth      0.014337   0.056920   0.252  0.80125    
## regionwest       0.080490   0.058784   1.369  0.17156    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4607 on 478 degrees of freedom
## Multiple R-squared:  0.2721, Adjusted R-squared:  0.2523 
## F-statistic: 13.75 on 13 and 478 DF,  p-value: < 2.2e-16
drop_vars_old <- names(Filter(function(x) length(unique(x)) <= 1, old[categorical_vars]))
old_clean <- old[, !(names(old) %in% drop_vars_old)]

model_old_full <- lm(log(wage) ~ . - wage - age, data = old_clean)
## Warning in terms.formula(formula, data = data): 'varlist' has changed (from
## nvar=11) to new 12 after EncodeVars() -- should no longer happen!
summary(model_old_full)
## 
## Call:
## lm(formula = log(wage) ~ . - wage - age, data = old_clean)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.69157 -0.30955 -0.00773  0.30861  1.30329 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      2.112012   0.158035  13.364  < 2e-16 ***
## educ             0.153919   0.012980  11.858  < 2e-16 ***
## hrswork          0.003061   0.002275   1.346 0.178857    
## insureYes        0.207256   0.058639   3.534 0.000437 ***
## metro1           0.133909   0.051927   2.579 0.010125 *  
## nchild          -0.014364   0.018738  -0.767 0.443580    
## union1           0.052524   0.053371   0.984 0.325408    
## raceBlack       -0.022868   0.112400  -0.203 0.838843    
## raceWhite        0.087461   0.093605   0.934 0.350448    
## marital1         0.139317   0.059422   2.345 0.019340 *  
## marital2         0.149507   0.071389   2.094 0.036610 *  
## regionnortheast  0.052852   0.058835   0.898 0.369341    
## regionsouth      0.055216   0.050949   1.084 0.278866    
## regionwest       0.172186   0.057079   3.017 0.002652 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4975 on 675 degrees of freedom
## Multiple R-squared:  0.2613, Adjusted R-squared:  0.2471 
## F-statistic: 18.37 on 13 and 675 DF,  p-value: < 2.2e-16

Multiple linear regression models were constructed separately for the young (age < 40) and old (age ≥ 40) subsets, with the natural logarithm of wage as the dependent variable and all other relevant variables (excluding wage and age) as predictors (Zeini, Okasha and Soliman, 2023). Categorical variables such as gender, race, marital status, region, metropolitan area, union membership, and insurance status were appropriately converted into factors to ensure correct model interpretation. Comparing the two models, differences in coefficient estimates were observed, reflecting how various factors influence wages differently across age groups. The R² values of the multiple regression models were higher than those from the simple linear regressions, indicating improved explanatory power when incorporating multiple predictors. However, the increase in R² also raises the possibility of overfitting, especially with many predictors relative to sample size. Reduced models with fewer variables are often preferable to the full models due to several reasons. Firstly, multicollinearity between predictors can distort coefficient estimates and reduce model reliability. Secondly, simpler models enhance interpretability, making it easier to understand the key factors affecting wages. Lastly, parsimony reduces the risk of overfitting and improves the model’s generalizability to other data sets. Therefore, variable selection and model refinement are important steps in developing effective regression models.

Conclusion

This analysis revealed notable wage disparities influenced by factors such as age, gender, and region, with older individuals generally earning higher wages. Insurance coverage showed some variation by gender, and family size impacted wage estimates. The regression models, both simple and multiple, effectively captured key relationships, with multiple regression improving explanatory power. However, limitations included small sample sizes for some subgroups and potential multicollinearity among predictors, which may affect model reliability. Future work could focus on refining models and incorporating additional variables to better understand wage determinants.

References

Bakhtsiyarava, M., Williams, T.G., Verdin, A. and Guikema, S.D., 2021. A nonparametric analysis of household-level food insecurity and its determinant factors: exploratory study in Ethiopia and Nigeria. Food Security, 13, pp.55-70.

José, T.R., Jhoset, Y.A., Soria, J.J. and Saboya, N., 2024, April. Machine Learning Models for Salary Prediction in Peruvian Teachers of Regular Basic Education. In Computer Science On-line Conference (pp. 534-552). Cham: Springer Nature Switzerland.

Kamruzzaman, M., Horowitz, M., Rahman, M.S., Deshmukh, H., Jones, K.L. and Marathe, C.S., 2025. Glycemic control is worse in rural compared to urban type 2 diabetes in Bangladesh, irrespective of food security status. Journal of Diabetes Investigation.

Nicodemo, C. and Satorra, A., 2022. Exploratory data analysis on large data sets: The example of salary variation in Spanish Social Security Data. BRQ Business Research Quarterly, 25(3), pp.283-294.

Tennhardt, L.M., Lambin, E.F., Curran, M. and Schader, C., 2023. Implementation of sustainable farming practices by cocoa farmers in Ecuador and Uganda: the influence of value chain factors. Frontiers in Sustainable Food Systems, 7, p.1167683.

Zeini, N.T., Okasha, A.E. and Soliman, A.S., 2023. Exploring and measuring quality of life determinants of wage workers in Egypt: a structural equation modelling approach. Social Indicators Research, 170(2), pp.339-374.