R Notebook

#list variable

colnames(Fulldat_Pheno)

 [1] "SEQN"                                  "chronological_age"                     "Gender"                                "Race"                                 
 [5] "Pregnancy"                             "Marital_Status"                        "Ratio_income_poverty"                  "Interview_Weight"                     
 [9] "MEC_Weight"                            "psu"                                   "Strata"                                "BMI"                                  
[13] "Vitamin_A"                             "Vitamin_C"                             "Vitamin_E"                             "Zinc"                                 
[17] "Selenium"                              "sleep_disorders"                       "Smoked_days"                           "now_smoke"                            
[21] "quit_smoking"                          "Avg_alcohol_drinks"                    "equipment_walk"                        "walk_difficulty"                      
[25] "had_cancer"                            "weight_2"                              "Perfluorohexane_sulfonic_acid"         "Perfluorohexane_sulfonic_acid_comment"
[29] "Perfluorononanoic_acid"                "Perfluorononanoic_acid_comment"        "perfluorooctanoic_acid"                "perfluorooctanoic_acid_comment"       
[33] "perfluorooctane_sulfonic_acid"         "perfluorooctane_sulfonic_acid_comment" "White_blood_cell_count"                "Lymphocyte_percent"                   
[37] "Red_cell_distribution_width"           "Mean_cell_volume"                      "Albumin"                               "Creatinine"                           
[41] "Glucose_serum"                         "Alkaline_phosphotase"                  "xb"                                    "Phenotypic_Age"                       
[45] "cate_age"                              "age_binary"                           

#Examine the pfas and Phenotypic_Age

ggplot(Fulldat_Pheno, aes(x = Perfluorohexane_sulfonic_acid)) +
  geom_histogram(binwidth = 100, color = "skyblue", fill = "red", alpha = 0.7) +
  labs(title = "Distribution of Perfluorohexane_sulfonic_acid",
       x = "Perfluorohexane_sulfonic_acid",
       y = "Frequency") +
  theme_minimal()

ggplot(Fulldat_Pheno, aes(x = Perfluorohexane_sulfonic_acid, y = Phenotypic_Age)) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE) + 
  labs(x = "Perfluorohexane_sulfonic_acid", y = "Phenotypic_Age", title = "Scatter Plot Perfluorohexane_sulfonic_acid vs Phenotypic_Age with Regression Line")

# Define breaks for age groups
breaks <- c(20, 45, 65, Inf)  # Breaks represent the age boundaries

# Define labels for the age groups
labels <- c("1", "2", "3")

# Categorize chronological_age into groups and assign custom labels
Fulldat_Pheno$cate_age <- cut(Fulldat_Pheno$chronological_age, breaks = breaks, labels = labels, include.lowest = TRUE)

#sample density curves of pfas concentrations among accelerated and delayed age

Fulldat_Pheno <- Fulldat_Pheno %>% mutate(age_binary = case_when(
  Phenotypic_Age-chronological_age >= 0 ~"accelerated",
  Phenotypic_Age-chronological_age < 0 ~"delayed"
))

library(ggplot2)
cols <- c("#F76D5E", "#72D8FF")

#Perfluorohexane_sulfonic_acid
ggplot(Fulldat_Pheno, aes(x = Perfluorohexane_sulfonic_acid, colour = age_binary)) +
  geom_density(lwd = 1.2, linetype = 1) + 
  scale_color_manual(values = cols)

#Perfluorohexane_sulfonic_acid_comment
ggplot(Fulldat_Pheno, aes(x = Perfluorononanoic_acid, colour = age_binary)) +
  geom_density(lwd = 1.2, linetype = 1) + 
  scale_color_manual(values = cols)

#Perfluorononanoic_acid
ggplot(Fulldat_Pheno, aes(x = perfluorooctanoic_acid, colour = age_binary)) +
  geom_density(lwd = 1.2, linetype = 1) + 
  scale_color_manual(values = cols)

#Perfluorononanoic_acid_comment
ggplot(Fulldat_Pheno, aes(x = perfluorooctane_sulfonic_acid, colour = age_binary)) +
  geom_density(lwd = 1.2, linetype = 1) + 
  scale_color_manual(values = cols)

#Main model of regression for association, and adjust for covariates (Table 2) #Perfluorohexane_sulfonic_acid

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid_comment, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            45.4975     0.2657 171.230   <2e-16 ***
Perfluorohexane_sulfonic_acid_comment   5.2534     2.1305   2.466   0.0149 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 377.0007)

Number of Fisher Scoring iterations: 2

                                          2.5 %    97.5 %
(Intercept)                           44.971890 46.023018
Perfluorohexane_sulfonic_acid_comment  1.039329  9.467451

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            50.2773     1.4506  34.660  < 2e-16 ***
Perfluorohexane_sulfonic_acid_comment   5.9817     2.2028   2.716  0.00752 ** 
Gender                                 -0.3176     0.4117  -0.771  0.44188    
Race                                    1.2279     0.1857   6.613 9.05e-10 ***
Marital_Status                         -2.9136     0.3698  -7.879 1.20e-12 ***
Ratio_income_poverty                   -0.3267     0.1731  -1.887  0.06146 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 336.8121)

Number of Fisher Scoring iterations: 2

                                           2.5 %      97.5 %
(Intercept)                           47.4072877 53.14741148
Perfluorohexane_sulfonic_acid_comment  1.6234305 10.33995224
Gender                                -1.1321702  0.49698242
Race                                   0.8605062  1.59519676
Marital_Status                        -3.6452647 -2.18194568
Ratio_income_poverty                  -0.6692319  0.01591455

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                        Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            7.853e+01  5.157e+00  15.227  < 2e-16 ***
Perfluorohexane_sulfonic_acid_comment  1.475e+01  3.398e+00   4.342 3.56e-05 ***
Gender                                -9.155e-01  9.737e-01  -0.940  0.34950    
Race                                   1.482e+00  4.314e-01   3.436  0.00088 ***
Marital_Status                        -3.117e+00  3.603e-01  -8.650 1.36e-13 ***
Ratio_income_poverty                  -1.089e+00  3.932e-01  -2.770  0.00675 ** 
BMI                                    2.014e-01  7.634e-02   2.638  0.00977 ** 
sleep_disorders                       -2.740e+00  1.133e+00  -2.418  0.01755 *  
quit_smoking                           3.096e-04  2.797e-05  11.069  < 2e-16 ***
Avg_alcohol_drinks                     5.775e+00  1.210e+00   4.772 6.67e-06 ***
had_cancer                            -1.480e+01  1.597e+00  -9.266 6.68e-15 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 305.0519)

Number of Fisher Scoring iterations: 2

                                              2.5 %        97.5 %
(Intercept)                            6.828635e+01  8.876488e+01
Perfluorohexane_sulfonic_acid_comment  8.005694e+00  2.149824e+01
Gender                                -2.848896e+00  1.017804e+00
Race                                   6.258390e-01  2.338989e+00
Marital_Status                        -3.832281e+00 -2.401372e+00
Ratio_income_poverty                  -1.869679e+00 -3.084527e-01
BMI                                    4.978659e-02  3.529389e-01
sleep_disorders                       -4.989964e+00 -4.897551e-01
quit_smoking                           2.540562e-04  3.651262e-04
Avg_alcohol_drinks                     3.371952e+00  8.177537e+00
had_cancer                            -1.796713e+01 -1.162605e+01

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   45.06187    0.32203 139.932  < 2e-16 ***
Perfluorohexane_sulfonic_acid  0.22852    0.08424   2.713  0.00756 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 376.9572)

Number of Fisher Scoring iterations: 2

                                    2.5 %     97.5 %
(Intercept)                   44.42491129 45.6988252
Perfluorohexane_sulfonic_acid  0.06189921  0.3951471

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                   49.56194    1.46138  33.914  < 2e-16 ***
Perfluorohexane_sulfonic_acid  0.26550    0.08535   3.111   0.0023 ** 
Gender                        -0.04536    0.42646  -0.106   0.9155    
Race                           1.20684    0.18319   6.588 1.03e-09 ***
Marital_Status                -2.92097    0.37017  -7.891 1.12e-12 ***
Ratio_income_poverty          -0.36969    0.17339  -2.132   0.0349 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 336.7355)

Number of Fisher Scoring iterations: 2

                                    2.5 %      97.5 %
(Intercept)                   46.67056007 52.45332776
Perfluorohexane_sulfonic_acid  0.09662619  0.43436501
Gender                        -0.88913426  0.79840851
Race                           0.84438890  1.56929608
Marital_Status                -3.65336306 -2.18858112
Ratio_income_poverty          -0.71273810 -0.02663413

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorohexane_sulfonic_acid + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                Estimate Std. Error t value Pr(>|t|)    
(Intercept)                    7.672e+01  5.511e+00  13.921  < 2e-16 ***
Perfluorohexane_sulfonic_acid  5.127e-01  2.134e-01   2.402  0.01826 *  
Gender                        -3.950e-01  1.023e+00  -0.386  0.70042    
Race                           1.361e+00  4.276e-01   3.184  0.00197 ** 
Marital_Status                -3.094e+00  3.601e-01  -8.594 1.78e-13 ***
Ratio_income_poverty          -1.112e+00  3.857e-01  -2.883  0.00488 ** 
BMI                            2.144e-01  7.670e-02   2.796  0.00628 ** 
sleep_disorders               -2.848e+00  1.121e+00  -2.541  0.01270 *  
quit_smoking                   3.060e-04  2.723e-05  11.240  < 2e-16 ***
Avg_alcohol_drinks             5.955e+00  1.226e+00   4.856 4.77e-06 ***
had_cancer                    -1.479e+01  1.645e+00  -8.987 2.62e-14 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 304.981)

Number of Fisher Scoring iterations: 2

                                      2.5 %        97.5 %
(Intercept)                    6.577500e+01  8.765862e+01
Perfluorohexane_sulfonic_acid  8.892741e-02  9.365442e-01
Gender                        -2.427095e+00  1.637121e+00
Race                           5.124391e-01  2.210311e+00
Marital_Status                -3.809366e+00 -2.379554e+00
Ratio_income_poverty          -1.877649e+00 -3.461801e-01
BMI                            6.212809e-02  3.666967e-01
sleep_disorders               -5.072929e+00 -6.223237e-01
quit_smoking                   2.519642e-04  3.600823e-04
Avg_alcohol_drinks             3.519771e+00  8.389605e+00
had_cancer                    -1.805244e+01 -1.151897e+01

#“Perfluorononanoic_acid” “Perfluorononanoic_acid_comment”

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid_comment, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     45.5414     0.2615 174.168   <2e-16 ***
Perfluorononanoic_acid_comment   1.2082     2.5944   0.466    0.642    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 377.3081)

Number of Fisher Scoring iterations: 2

                                   2.5 %    97.5 %
(Intercept)                    45.024247 46.058640
Perfluorononanoic_acid_comment -3.923447  6.339837

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     50.4123     1.4495  34.780  < 2e-16 ***
Perfluorononanoic_acid_comment  -0.1223     2.5110  -0.049   0.9612    
Gender                          -0.2958     0.4101  -0.721   0.4720    
Race                             1.2229     0.1839   6.650 7.54e-10 ***
Marital_Status                  -2.9142     0.3698  -7.881 1.18e-12 ***
Ratio_income_poverty            -0.3530     0.1733  -2.037   0.0437 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 337.238)

Number of Fisher Scoring iterations: 2

                                    2.5 %      97.5 %
(Intercept)                    47.5444365 53.28010557
Perfluorononanoic_acid_comment -5.0903542  4.84577011
Gender                         -1.1072206  0.51559131
Race                            0.8590130  1.58672663
Marital_Status                 -3.6458461 -2.18263403
Ratio_income_poverty           -0.6959141 -0.01015699

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     7.815e+01  5.256e+00  14.868  < 2e-16 ***
Perfluorononanoic_acid_comment -7.871e+00  4.689e+00  -1.678  0.09660 .  
Gender                         -7.635e-01  9.504e-01  -0.803  0.42383    
Race                            1.416e+00  4.272e-01   3.314  0.00131 ** 
Marital_Status                 -3.169e+00  3.543e-01  -8.946 3.20e-14 ***
Ratio_income_poverty           -1.157e+00  3.943e-01  -2.934  0.00421 ** 
BMI                             2.148e-01  7.719e-02   2.783  0.00651 ** 
sleep_disorders                -2.678e+00  1.136e+00  -2.358  0.02045 *  
quit_smoking                    3.112e-04  2.895e-05  10.750  < 2e-16 ***
Avg_alcohol_drinks              5.964e+00  1.248e+00   4.778 6.52e-06 ***
had_cancer                     -1.468e+01  1.606e+00  -9.139 1.25e-14 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 305.8616)

Number of Fisher Scoring iterations: 2

                                       2.5 %        97.5 %
(Intercept)                     67.713073728  8.858582e+01
Perfluorononanoic_acid_comment -17.181591762  1.440350e+00
Gender                          -2.650618824  1.123632e+00
Race                             0.567440494  2.264067e+00
Marital_Status                  -3.872912412 -2.466056e+00
Ratio_income_poverty            -1.939889398 -3.739266e-01
BMI                              0.061579130  3.681027e-01
sleep_disorders                 -4.932421575 -4.229425e-01
quit_smoking                     0.000253748  3.687221e-04
Avg_alcohol_drinks               3.485471884  8.442395e+00
had_cancer                     -17.864999070 -1.148764e+01

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                       Estimate Std. Error t value Pr(>|t|)    
(Intercept)             44.6374     0.3784 117.977  < 2e-16 ***
Perfluorononanoic_acid   0.8370     0.2218   3.774 0.000241 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 376.423)

Number of Fisher Scoring iterations: 2

                            2.5 %    97.5 %
(Intercept)            43.8890251 45.385773
Perfluorononanoic_acid  0.3983308  1.275581

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid + Gender + 
    Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                       Estimate Std. Error t value Pr(>|t|)    
(Intercept)             49.7547     1.4551  34.192  < 2e-16 ***
Perfluorononanoic_acid   0.7136     0.2264   3.152  0.00202 ** 
Gender                  -0.1994     0.4160  -0.479  0.63249    
Race                     1.1582     0.1869   6.198 7.14e-09 ***
Marital_Status          -2.9157     0.3702  -7.877 1.21e-12 ***
Ratio_income_poverty    -0.3803     0.1746  -2.178  0.03119 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 336.5506)

Number of Fisher Scoring iterations: 2

                            2.5 %      97.5 %
(Intercept)            46.8757035 52.63376825
Perfluorononanoic_acid  0.2656063  1.16153040
Gender                 -1.0224187  0.62360929
Race                    0.7884467  1.52787727
Marital_Status         -3.6481541 -2.18333322
Ratio_income_poverty   -0.7256677 -0.03490319

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid + Gender + 
    Race + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + 
    Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + 
    had_cancer, design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                         Estimate Std. Error t value Pr(>|t|)    
(Intercept)             7.782e+01  5.338e+00  14.579  < 2e-16 ***
Perfluorononanoic_acid  6.655e-01  5.132e-01   1.297  0.19784    
Gender                 -7.893e-01  9.966e-01  -0.792  0.43036    
Race                    1.331e+00  4.370e-01   3.046  0.00301 ** 
Marital_Status         -3.105e+00  3.614e-01  -8.592 1.81e-13 ***
Ratio_income_poverty   -1.108e+00  3.931e-01  -2.818  0.00589 ** 
BMI                     2.109e-01  7.612e-02   2.770  0.00675 ** 
sleep_disorders        -2.767e+00  1.126e+00  -2.458  0.01582 *  
quit_smoking            3.107e-04  2.799e-05  11.099  < 2e-16 ***
Avg_alcohol_drinks      5.889e+00  1.212e+00   4.858 4.72e-06 ***
had_cancer             -1.476e+01  1.625e+00  -9.084 1.63e-14 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 306.5129)

Number of Fisher Scoring iterations: 2

                               2.5 %        97.5 %
(Intercept)             6.722415e+01  8.842209e+01
Perfluorononanoic_acid -3.533787e-01  1.684429e+00
Gender                 -2.768086e+00  1.189481e+00
Race                    4.634000e-01  2.198819e+00
Marital_Status         -3.822612e+00 -2.387453e+00
Ratio_income_poverty   -1.888041e+00 -3.271788e-01
BMI                     5.975077e-02  3.620382e-01
sleep_disorders        -5.001788e+00 -5.314698e-01
quit_smoking            2.550832e-04  3.662334e-04
Avg_alcohol_drinks      3.482350e+00  8.295721e+00
had_cancer             -1.798581e+01 -1.153376e+01

#“perfluorooctanoic_acid” “perfluorooctanoic_acid_comment”

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctanoic_acid_comment, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     43.4970     0.3398 127.989  < 2e-16 ***
perfluorooctanoic_acid_comment   4.8208     0.6522   7.391 2.29e-11 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 368.6468)

Number of Fisher Scoring iterations: 2

                                   2.5 %    97.5 %
(Intercept)                    42.824027 44.170014
perfluorooctanoic_acid_comment  3.529163  6.112411

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctanoic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     48.5424     1.5971  30.393  < 2e-16 ***
perfluorooctanoic_acid_comment   3.1479     0.6242   5.043 1.74e-06 ***
Gender                          -0.3302     0.4575  -0.722   0.4719    
Race                             1.2890     0.2016   6.393 3.70e-09 ***
Marital_Status                  -2.8443     0.3903  -7.287 4.45e-11 ***
Ratio_income_poverty            -0.3586     0.1767  -2.030   0.0447 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 330.7506)

Number of Fisher Scoring iterations: 2

                                    2.5 %       97.5 %
(Intercept)                    45.3784365 51.706315821
perfluorooctanoic_acid_comment  1.9114453  4.384409757
Gender                         -1.2365070  0.576118449
Race                            0.8895389  1.688413922
Marital_Status                 -3.6175446 -2.071017131
Ratio_income_poverty           -0.7085828 -0.008648249

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctanoic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)                     7.819e+01  6.022e+00  12.985  < 2e-16 ***
perfluorooctanoic_acid_comment  1.786e+00  1.205e+00   1.482 0.142284    
Gender                         -1.031e+00  1.059e+00  -0.973 0.333315    
Race                            1.345e+00  4.793e-01   2.806 0.006318 ** 
Marital_Status                 -3.293e+00  3.767e-01  -8.741 3.12e-13 ***
Ratio_income_poverty           -1.074e+00  3.979e-01  -2.700 0.008491 ** 
BMI                             1.946e-01  8.413e-02   2.314 0.023286 *  
sleep_disorders                -2.881e+00  1.244e+00  -2.316 0.023140 *  
quit_smoking                    2.891e-04  3.467e-05   8.337 1.92e-12 ***
Avg_alcohol_drinks              5.953e+00  1.539e+00   3.868 0.000224 ***
had_cancer                     -1.475e+01  1.816e+00  -8.126 4.95e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 306.8411)

Number of Fisher Scoring iterations: 2

                                       2.5 %        97.5 %
(Intercept)                     6.620563e+01  9.017810e+01
perfluorooctanoic_acid_comment -6.125094e-01  4.184327e+00
Gender                         -3.137935e+00  1.076757e+00
Race                            3.907639e-01  2.298814e+00
Marital_Status                 -4.042829e+00 -2.543076e+00
Ratio_income_poverty           -1.866186e+00 -2.821458e-01
BMI                             2.719110e-02  3.620871e-01
sleep_disorders                -5.356931e+00 -4.052082e-01
quit_smoking                    2.200367e-04  3.580638e-04
Avg_alcohol_drinks              2.889848e+00  9.015653e+00
had_cancer                     -1.836698e+01 -1.113947e+01

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                       Estimate Std. Error t value Pr(>|t|)    
(Intercept)             44.6374     0.3784 117.977  < 2e-16 ***
Perfluorononanoic_acid   0.8370     0.2218   3.774 0.000241 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 376.423)

Number of Fisher Scoring iterations: 2

                            2.5 %    97.5 %
(Intercept)            43.8890251 45.385773
Perfluorononanoic_acid  0.3983308  1.275581

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid + Gender + 
    Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                       Estimate Std. Error t value Pr(>|t|)    
(Intercept)             49.7547     1.4551  34.192  < 2e-16 ***
Perfluorononanoic_acid   0.7136     0.2264   3.152  0.00202 ** 
Gender                  -0.1994     0.4160  -0.479  0.63249    
Race                     1.1582     0.1869   6.198 7.14e-09 ***
Marital_Status          -2.9157     0.3702  -7.877 1.21e-12 ***
Ratio_income_poverty    -0.3803     0.1746  -2.178  0.03119 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 336.5506)

Number of Fisher Scoring iterations: 2

                            2.5 %      97.5 %
(Intercept)            46.8757035 52.63376825
Perfluorononanoic_acid  0.2656063  1.16153040
Gender                 -1.0224187  0.62360929
Race                    0.7884467  1.52787727
Marital_Status         -3.6481541 -2.18333322
Ratio_income_poverty   -0.7256677 -0.03490319

Call:
svyglm(formula = Phenotypic_Age ~ Perfluorononanoic_acid + Gender + 
    Race + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + 
    Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + 
    had_cancer, design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                         Estimate Std. Error t value Pr(>|t|)    
(Intercept)             7.782e+01  5.338e+00  14.579  < 2e-16 ***
Perfluorononanoic_acid  6.655e-01  5.132e-01   1.297  0.19784    
Gender                 -7.893e-01  9.966e-01  -0.792  0.43036    
Race                    1.331e+00  4.370e-01   3.046  0.00301 ** 
Marital_Status         -3.105e+00  3.614e-01  -8.592 1.81e-13 ***
Ratio_income_poverty   -1.108e+00  3.931e-01  -2.818  0.00589 ** 
BMI                     2.109e-01  7.612e-02   2.770  0.00675 ** 
sleep_disorders        -2.767e+00  1.126e+00  -2.458  0.01582 *  
quit_smoking            3.107e-04  2.799e-05  11.099  < 2e-16 ***
Avg_alcohol_drinks      5.889e+00  1.212e+00   4.858 4.72e-06 ***
had_cancer             -1.476e+01  1.625e+00  -9.084 1.63e-14 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 306.5129)

Number of Fisher Scoring iterations: 2

                               2.5 %        97.5 %
(Intercept)             6.722415e+01  8.842209e+01
Perfluorononanoic_acid -3.533787e-01  1.684429e+00
Gender                 -2.768086e+00  1.189481e+00
Race                    4.634000e-01  2.198819e+00
Marital_Status         -3.822612e+00 -2.387453e+00
Ratio_income_poverty   -1.888041e+00 -3.271788e-01
BMI                     5.975077e-02  3.620382e-01
sleep_disorders        -5.001788e+00 -5.314698e-01
quit_smoking            2.550832e-04  3.662334e-04
Avg_alcohol_drinks      3.482350e+00  8.295721e+00
had_cancer             -1.798581e+01 -1.153376e+01

#“perfluorooctane_sulfonic_acid” “perfluorooctane_sulfonic_acid_comment”

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid_comment, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            45.1300     0.2921 154.485   <2e-16 ***
perfluorooctane_sulfonic_acid_comment   8.2645     4.1835   1.975   0.0505 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 373.6055)

Number of Fisher Scoring iterations: 2

                                            2.5 %   97.5 %
(Intercept)                           44.55146084 45.70846
perfluorooctane_sulfonic_acid_comment -0.02003231 16.54910

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            49.6362     1.6016  30.991  < 2e-16 ***
perfluorooctane_sulfonic_acid_comment   9.6344     5.2158   1.847   0.0673 .  
Gender                                 -0.3451     0.4578  -0.754   0.4525    
Race                                    1.3338     0.1969   6.773 5.80e-10 ***
Marital_Status                         -2.9403     0.3994  -7.362 3.03e-11 ***
Ratio_income_poverty                   -0.3481     0.1821  -1.911   0.0585 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 332.5743)

Number of Fisher Scoring iterations: 2

                                           2.5 %      97.5 %
(Intercept)                           46.4633496 52.80905212
perfluorooctane_sulfonic_acid_comment -0.6980838 19.96687482
Gender                                -1.2519299  0.56172801
Race                                   0.9437016  1.72397764
Marital_Status                        -3.7314691 -2.14916351
Ratio_income_poverty                  -0.7088672  0.01268752

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid_comment + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                        Estimate Std. Error t value Pr(>|t|)    
(Intercept)                            7.928e+01  5.907e+00  13.421  < 2e-16 ***
perfluorooctane_sulfonic_acid_comment  1.955e+00  8.065e+00   0.242 0.809060    
Gender                                -9.339e-01  1.079e+00  -0.865 0.389401    
Race                                   1.390e+00  4.734e-01   2.936 0.004350 ** 
Marital_Status                        -3.376e+00  4.007e-01  -8.426 1.29e-12 ***
Ratio_income_poverty                  -1.080e+00  4.013e-01  -2.692 0.008669 ** 
BMI                                    2.052e-01  8.167e-02   2.512 0.014037 *  
sleep_disorders                       -3.008e+00  1.246e+00  -2.414 0.018112 *  
quit_smoking                           3.012e-04  3.324e-05   9.063 7.34e-14 ***
Avg_alcohol_drinks                     5.591e+00  1.486e+00   3.762 0.000322 ***
had_cancer                            -1.481e+01  1.831e+00  -8.090 5.83e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 307.6054)

Number of Fisher Scoring iterations: 2

                                              2.5 %        97.5 %
(Intercept)                            6.752297e+01  9.103943e+01
perfluorooctane_sulfonic_acid_comment -1.409706e+01  1.800763e+01
Gender                                -3.081612e+00  1.213880e+00
Race                                   4.476930e-01  2.332064e+00
Marital_Status                        -4.173432e+00 -2.578424e+00
Ratio_income_poverty                  -1.878813e+00 -2.814733e-01
BMI                                    4.260483e-02  3.677075e-01
sleep_disorders                       -5.487817e+00 -5.272295e-01
quit_smoking                           2.350757e-04  3.673987e-04
Avg_alcohol_drinks                     2.632925e+00  8.549037e+00
had_cancer                            -1.845341e+01 -1.116557e+01

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid, 
    design = des, family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                    43.0812     0.3926 109.739  < 2e-16 ***
perfluorooctane_sulfonic_acid   0.1562     0.0196   7.971 1.12e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 368.2211)

Number of Fisher Scoring iterations: 2

                                   2.5 %     97.5 %
(Intercept)                   42.3037475 43.8585774
perfluorooctane_sulfonic_acid  0.1173858  0.1949941

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid + 
    Gender + Race + Marital_Status + Ratio_income_poverty, design = des, 
    family = "gaussian", data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients:
                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                    47.5011     1.6129  29.451  < 2e-16 ***
perfluorooctane_sulfonic_acid   0.1322     0.0192   6.885 3.34e-10 ***
Gender                          0.3089     0.4690   0.659   0.5115    
Race                            1.1553     0.1969   5.866 4.44e-08 ***
Marital_Status                 -2.9317     0.3985  -7.357 3.11e-11 ***
Ratio_income_poverty           -0.3780     0.1806  -2.092   0.0386 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 328.8112)

Number of Fisher Scoring iterations: 2

                                    2.5 %      97.5 %
(Intercept)                   44.30600893 50.69617633
perfluorooctane_sulfonic_acid  0.09416215  0.17024027
Gender                        -0.62015706  1.23788344
Race                           0.76516030  1.54538803
Marital_Status                -3.72110932 -2.14234470
Ratio_income_poverty          -0.73575239 -0.02013906

Call:
svyglm(formula = Phenotypic_Age ~ perfluorooctane_sulfonic_acid + 
    Gender + Race + Marital_Status + Ratio_income_poverty + BMI + 
    sleep_disorders + Smoked_days + now_smoke + quit_smoking + 
    Avg_alcohol_drinks + had_cancer, design = des, family = "gaussian", 
    data = Fulldat_Pheno)

Survey design:
svydesign(id = ~psu, strata = ~Strata, weights = ~weight_2, nest = TRUE, 
    data = Fulldat_Pheno)

Coefficients: (2 not defined because of singularities)
                                Estimate Std. Error t value Pr(>|t|)    
(Intercept)                    7.642e+01  6.089e+00  12.550  < 2e-16 ***
perfluorooctane_sulfonic_acid  1.821e-01  3.734e-02   4.876 5.49e-06 ***
Gender                         5.954e-02  1.111e+00   0.054 0.957399    
Race                           1.002e+00  4.635e-01   2.162 0.033670 *  
Marital_Status                -3.297e+00  4.033e-01  -8.176 3.96e-12 ***
Ratio_income_poverty          -1.028e+00  3.879e-01  -2.650 0.009725 ** 
BMI                            2.250e-01  8.148e-02   2.762 0.007149 ** 
sleep_disorders               -3.069e+00  1.243e+00  -2.468 0.015731 *  
quit_smoking                   2.995e-04  3.168e-05   9.457 1.25e-14 ***
Avg_alcohol_drinks             5.153e+00  1.405e+00   3.667 0.000444 ***
had_cancer                    -1.464e+01  1.827e+00  -8.015 8.17e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 301.645)

Number of Fisher Scoring iterations: 2

                                      2.5 %        97.5 %
(Intercept)                    6.430133e+01  88.542926525
perfluorooctane_sulfonic_acid  1.077436e-01   0.256374471
Gender                        -2.151987e+00   2.271065687
Race                           7.936354e-02   1.924426690
Marital_Status                -4.100053e+00  -2.494555627
Ratio_income_poverty          -1.800105e+00  -0.255781772
BMI                            6.284244e-02   0.387221415
sleep_disorders               -5.543968e+00  -0.594275299
quit_smoking                   2.365019e-04   0.000362604
Avg_alcohol_drinks             2.355575e+00   7.950191172
had_cancer                    -1.828070e+01 -11.007005084

run cubic spline model for non-linear regression(Figure 2)

#library("rcssci")
#rcssci_linear(data = Fulldat_Pheno, y = "Phenotypic_Age", x = "Perfluorohexane_sulfonic_acid", covs=c("Gender", "Race", "BMI","had_cancer"), prob = 0.1,ref.zero=FALSE,
#              filepath = "C:/Users/HKUSCM/Documents")
#

#check effect modifiers (subgroup Figure 3) #```{r echo=FALSE,message=FALSE,warning=TRUE} #Subgroup analysis for gender, race, BMI, income, smoking and cancer #for subgroup analysis by gender men <- Fulldat_Pheno[Fulldat_Pheno\(Gender == 1, ] women <- Fulldat_Pheno[Fulldat_Pheno\)Gender == 2, ] desmen <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = men) #select men deswomen <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = women) #select women

model_sex <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*Gender, data = Fulldat_Pheno) summary(model_sex)

model_men <- svyglm(Phenotypic_Age ~ categorize_MET, data = men, design = desmen, family = “gaussian”) summary(model_men) confint(model_men)[“categorize_MET2”, ]
model_women <- svyglm(Phenotypic_Age ~ categorize_MET, data = women, design = deswomen, family = “gaussian”) summary(model_women) confint(model_women)[“categorize_MET2”, ]

#for subgroup analysis by Race Mexican <- Fulldat_Pheno[Fulldat_Pheno\(Race == 1, ] Other_Hispanic <- Fulldat_Pheno[Fulldat_Pheno\)Race == 2, ] Non_Hispanic_white <- Fulldat_Pheno[Fulldat_Pheno\(Race == 3, ] Non_Hispanic_Black <- Fulldat_Pheno[Fulldat_Pheno\)Race == 4, ] Other_Race <- Fulldat_Pheno[Fulldat_Pheno$Race == 5, ]

model_race <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*as.factor(Race), data = Fulldat_Pheno) summary(model_race)

des_Me <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Mexican) des_Hispanic <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Other_Hispanic) des_white <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Non_Hispanic_white) des_Black <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Non_Hispanic_Black) des_Other <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Other_Race)

model_Me <- svyglm(Phenotypic_Age ~ categorize_MET, data = Mexican, design = des_Me, family = “gaussian”) summary(model_Me) confint(model_Me)[“categorize_MET2”, ]
model_Hispanic <- svyglm(Phenotypic_Age ~ categorize_MET, data = Other_Hispanic, design = des_Hispanic, family = “gaussian”) summary(model_Hispanic) confint(model_Hispanic)[“categorize_MET2”, ]
model_white <- svyglm(Phenotypic_Age ~ categorize_MET, data = Non_Hispanic_white, design = des_white, family = “gaussian”) summary(model_white) confint(model_white)[“categorize_MET2”, ]
model_Black <- svyglm(Phenotypic_Age ~ categorize_MET, data = Non_Hispanic_Black, design = des_Black, family = “gaussian”) summary(model_Black) confint(model_Black)[“categorize_MET2”, ]
model_Other <- svyglm(Phenotypic_Age ~ categorize_MET, data = Other_Race, design = des_Other, family = “gaussian”) summary(model_Other) confint(model_Other)[“categorize_MET2”, ]

##for subgroup analysis by BMI Fulldat_Pheno\(BMI[is.na(Fulldat_Pheno\)BMI)] <- 0 BMI_1 <- Fulldat_Pheno[Fulldat_Pheno\(BMI < 25, ] BMI_2 <- Fulldat_Pheno[Fulldat_Pheno\)BMI >= 25 & Fulldat_Pheno\(BMI < 30, ] BMI_3 <- Fulldat_Pheno[Fulldat_Pheno\)BMI >= 30, ]

Fulldat_Pheno <- Fulldat_Pheno %>% mutate(BMI_cat = case_when( BMI >= 30 ~ “obesity”, BMI <= 25 ~ “normal”, TRUE ~ “overweight” ) ) model_BMI <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*BMI_cat, data = Fulldat_Pheno) summary(model_BMI)

desBMI_1 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = BMI_1) #select BMI < 25 desBMI_2 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = BMI_2) #select 25 =< BMI <30 desBMI_3 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = BMI_3) #select >= 30 model_BMI_1 <- svyglm(Phenotypic_Age ~ categorize_MET, data = BMI_1, design = desBMI_1, family = “gaussian”) summary(model_BMI_1) confint(model_BMI_1)[“categorize_MET2”, ]
model_BMI_2 <- svyglm(Phenotypic_Age ~ categorize_MET, data = BMI_2, design = desBMI_2, family = “gaussian”) summary(model_BMI_2) confint(model_BMI_2)[“categorize_MET2”, ]
model_BMI_3 <- svyglm(Phenotypic_Age ~ categorize_MET, data = BMI_3, design = desBMI_3, family = “gaussian”) summary(model_BMI_3) confint(model_BMI_3)[“categorize_MET2”, ]

##for subgroup analysis by income Fulldat_Pheno\(Ratio_income_poverty[is.na(Fulldat_Pheno\)Ratio_income_poverty)] <- 0 income_1 <- Fulldat_Pheno[Fulldat_Pheno\(Ratio_income_poverty <= 1, ] income_2 <- Fulldat_Pheno[Fulldat_Pheno\)Ratio_income_poverty > 1 & Fulldat_Pheno\(Ratio_income_poverty < 4, ] income_3 <- Fulldat_Pheno[Fulldat_Pheno\)Ratio_income_poverty >= 4, ]

Fulldat_Pheno <- Fulldat_Pheno %>% mutate(income_cat = case_when( Ratio_income_poverty >= 4 ~ “rich”, Ratio_income_poverty <= 1 ~ “poor”, TRUE ~ “average” ) )

model_income <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*as.factor(income_cat), data = Fulldat_Pheno) summary(model_income)

desincome_1 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = income_1) #select income <= 1 desincome_2 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = income_2) #select 1 < income < 4 desincome_3 <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = income_3) #select >= 4 model_income_1 <- svyglm(Phenotypic_Age ~ categorize_MET, data = income_1, design = desincome_1, family = “gaussian”) summary(model_income_1) confint(model_income_1)[“categorize_MET2”, ] model_income_2 <- svyglm(Phenotypic_Age ~ categorize_MET, data = income_2, design = desincome_2, family = “gaussian”) summary(model_income_2) confint(model_income_2)[“categorize_MET2”, ] model_income_3 <- svyglm(Phenotypic_Age ~ categorize_MET, data = income_3, design = desincome_3, family = “gaussian”) summary(model_income_3) confint(model_income_3)[“categorize_MET2”, ]

#for subgroup analysis by cancer cancer <- Fulldat_Pheno[Fulldat_Pheno\(had_cancer == 1, ] cancer\)had_cancer[is.na(cancer\(had_cancer)] <- 0 cancer <- cancer[cancer\)had_cancer != 0, , drop = FALSE]

non_cancer <- Fulldat_Pheno[Fulldat_Pheno\(had_cancer == 2, ] non_cancer\)had_cancer[is.na(non_cancer\(had_cancer)] <- 0 non_cancer <- non_cancer[non_cancer\)had_cancer != 0, , drop = FALSE]

model_cancer <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*had_cancer, data = Fulldat_Pheno) summary(model_cancer)

des_cancer <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = cancer) des_non_cancer <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = non_cancer)

model_cancer <- svyglm(Phenotypic_Age ~ categorize_MET, data = cancer, design = des_cancer, family = “gaussian”) summary(model_cancer) confint(model_cancer)[“categorize_MET2”, ] model_non_cancer <- svyglm(Phenotypic_Age ~ categorize_MET, data = non_cancer, design = des_non_cancer, family = “gaussian”) summary(model_non_cancer) confint(model_non_cancer)[“categorize_MET2”, ]

##for subgroup analysis by smoking Fulldat_Pheno\(quit_smoking[is.na(Fulldat_Pheno\)quit_smoking)] <- 0 Fulldat_Pheno\(now_smoke[is.na(Fulldat_Pheno\)now_smoke)] <- 0 current_smokers <- Fulldat_Pheno[Fulldat_Pheno\(Smoked_days == 1 & Fulldat_Pheno\)now_smoke == 1, ] former_smokers <- Fulldat_Pheno[Fulldat_Pheno\(Smoked_days == 1 & Fulldat_Pheno\)quit_smoking > 1, ] former_smokers\(psu[is.na(former_smokers\)psu)] <- 0 former_smokers <- former_smokers[former_smokers$psu != 0, , drop = FALSE]

model_smoke <- glm(Phenotypic_Age ~ categorize_MET + categorize_MET*now_smoke, data = Fulldat_Pheno) summary(model_smoke)

non_smokers <- Fulldat_Pheno[Fulldat_Pheno\(Smoked_days == 2, ] non_smokers\)Smoked_days[is.na(non_smokers\(Smoked_days)] <- 0 non_smokers <- non_smokers[non_smokers\)Smoked_days != 0, , drop = FALSE] #select hose who were considered current smokers smoked on a regular basis and had smoked at least 100 cigarettes in their lifetime. descurrent_smokers <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = current_smokers) #select Former smokers had smoked at least 100 cigarettes and had since quit. desformer_smokers <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = former_smokers) #select Non-smokers had either never smoked or smoked fewer than 100 cigarettes desnon_smokers <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = non_smokers)

model_current_smokers <- svyglm(Phenotypic_Age ~ categorize_MET, data = current_smokers, design = descurrent_smokers, family = “gaussian”) summary(model_current_smokers) confint(model_current_smokers)[“categorize_MET2”, ] model_former_smokers <- svyglm(Phenotypic_Age ~ categorize_MET, data = former_smokers, design = desformer_smokers, family = “gaussian”) summary(model_former_smokers) confint(model_former_smokers)[“categorize_MET2”, ] model_non_smokers <- svyglm(Phenotypic_Age ~ categorize_MET, data = non_smokers, design = desnon_smokers, family = “gaussian”) summary(model_non_smokers) confint(model_non_smokers)[“categorize_MET2”, ]

#sensitivity analysis without cancer patient only for MET (Supple table 1)
#```{r}
#subset cancer
cancer <- Fulldat_Pheno[Fulldat_Pheno$had_cancer == 1, ]
cancer$had_cancer[is.na(cancer$had_cancer)] <- 0
cancer <- cancer[cancer$had_cancer != 0, , drop = FALSE]

non_cancer <- Fulldat_Pheno[Fulldat_Pheno$had_cancer == 2, ]
non_cancer$had_cancer[is.na(non_cancer$had_cancer)] <- 0
non_cancer <- non_cancer[non_cancer$had_cancer != 0, , drop = FALSE]

des_cancer <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = cancer) 
des_non_cancer <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = non_cancer) 

model_cancer1 <- svyglm(Phenotypic_Age ~ Total_MET, data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer1)
confint(model_cancer1)

model_cancer2 <- svyglm(Phenotypic_Age ~ Total_MET + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                   , data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer2)
confint(model_cancer2)


model_cancer3 <- svyglm(Phenotypic_Age ~ Total_MET+ Gender + Race + CDAI + Marital_Status + Ratio_income_poverty 
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer, data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer3)
confint(model_cancer3)

model_cancer4 <- svyglm(Phenotypic_Age ~ categorize_MET, data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer4)
confint(model_cancer4)


model_cancer5 <- svyglm(Phenotypic_Age ~ categorize_MET + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                   , data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer5)
confint(model_cancer5)


model_cancer6 <- svyglm(Phenotypic_Age ~ categorize_MET+ Gender + Race + CDAI + Marital_Status + Ratio_income_poverty 
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer, data = non_cancer, design = des_non_cancer, family = "gaussian")
summary(model_cancer6)
confint(model_cancer6)

#sensitivity analysis for Main model of regression for association on PA, and adjust for covariates (Suppl Table 2) #```{r} Fulldat_Pheno\(moderate_total_week <- Fulldat_Pheno\)Days_moderate_recreational * Fulldat_Pheno\(Minutes_moderate_recreational Fulldat_Pheno\)vigorous_total_week <- Fulldat_Pheno\(Days_vigorous_recreational * Fulldat_Pheno\)Minutes_vigorous_recreational Fulldat_Pheno\(moderate_total_day <- Fulldat_Pheno\)Minutes_moderate_recreational Fulldat_Pheno\(vigorous_total_day <- Fulldat_Pheno\)Minutes_vigorous_recreational

Fulldat_Pheno\(moderate_cat_week <- ifelse(Fulldat_Pheno\)moderate_total_week > 150, “1”, “0”) Fulldat_Pheno\(vigorous_cat_week <- ifelse(Fulldat_Pheno\)vigorous_total_week > 75, “1”, “0”) Fulldat_Pheno\(moderate_cat_day <- ifelse(Fulldat_Pheno\)Minutes_moderate_recreational > 60, “1”, “0”) Fulldat_Pheno\(vigorous_cat_day <- ifelse(Fulldat_Pheno\)Minutes_vigorous_recreational > 30, “1”, “0”)

des <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Fulldat_Pheno)

#vigorous_total_week model_s1 <- svyglm(Phenotypic_Age ~ vigorous_total_week, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_s1) confint(model_s1) model_s2 <- svyglm(Phenotypic_Age ~ vigorous_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_s2) confint(model_s2) model_s3 <- svyglm(Phenotypic_Age ~ vigorous_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_s3) confint(model_s3)

#moderate_total_week model_ss1 <- svyglm(Phenotypic_Age ~ moderate_total_week, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ss1) confint(model_ss1) model_ss2 <- svyglm(Phenotypic_Age ~ moderate_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ss2) confint(model_ss2) model_ss3 <- svyglm(Phenotypic_Age ~ moderate_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ss3) confint(model_ss3)

#vigorous_total_day model_sss1 <- svyglm(Phenotypic_Age ~ vigorous_total_day, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sss1) confint(model_sss1) model_sss2 <- svyglm(Phenotypic_Age ~ vigorous_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sss2) confint(model_sss2) model_sss3 <- svyglm(Phenotypic_Age ~ vigorous_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sss3) confint(model_sss3)

#moderate_total_day model_ssss1 <- svyglm(Phenotypic_Age ~ moderate_total_day, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssss1) confint(model_ssss1) model_ssss2 <- svyglm(Phenotypic_Age ~ moderate_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssss2) confint(model_ssss2) model_ssss3 <- svyglm(Phenotypic_Age ~ moderate_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssss3) confint(model_ssss3)

#vigorous_cat_week model_sssss1 <- svyglm(Phenotypic_Age ~ vigorous_cat_week, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssss1) confint(model_sssss1) model_sssss2 <- svyglm(Phenotypic_Age ~ vigorous_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssss2) confint(model_sssss2) model_sssss3 <- svyglm(Phenotypic_Age ~ vigorous_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssss3) confint(model_sssss3)

#moderate_cat_week model_ssssss1 <- svyglm(Phenotypic_Age ~ moderate_cat_week, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssss1) confint(model_ssssss1) model_ssssss2 <- svyglm(Phenotypic_Age ~ moderate_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssss2) confint(model_ssssss2) model_ssssss3 <- svyglm(Phenotypic_Age ~ moderate_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssss3) confint(model_ssssss3)

#vigorous_cat_day model_sssssss1 <- svyglm(Phenotypic_Age ~ vigorous_cat_day, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssssss1) confint(model_sssssss1) model_sssssss2 <- svyglm(Phenotypic_Age ~ vigorous_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssssss2) confint(model_sssssss2) model_sssssss3 <- svyglm(Phenotypic_Age ~ vigorous_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_sssssss3) confint(model_sssssss3)

#moderate_cat_day model_ssssssss1 <- svyglm(Phenotypic_Age ~ moderate_cat_day, data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssssss1) confint(model_ssssssss1) model_ssssssss2 <- svyglm(Phenotypic_Age ~ moderate_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssssss2) confint(model_ssssssss2) model_ssssssss3 <- svyglm(Phenotypic_Age ~ moderate_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer , data = Fulldat_Pheno, design = des, family = “gaussian”) summary(model_ssssssss3) confint(model_ssssssss3)

#```{r}
des <- svydesign(id =~ psu, strata =~ Strata, weights =~ weight_2, nest = TRUE, data = Fulldat_Pheno)

#vigorous_total_week
model_s1 <- svyglm(Phenotypic_Age ~ vigorous_total_week, data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_s1)
confint(model_s1)
model_s2 <- svyglm(Phenotypic_Age ~ vigorous_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_s2)
confint(model_s2)
model_s3 <- svyglm(Phenotypic_Age ~ vigorous_total_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_s3)
confint(model_s3)


#vigorous_total_day
model_sss1 <- svyglm(Phenotypic_Age ~ vigorous_total_day, data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sss1)
confint(model_sss1)
model_sss2 <- svyglm(Phenotypic_Age ~ vigorous_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sss2)
confint(model_sss2)
model_sss3 <- svyglm(Phenotypic_Age ~ vigorous_total_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sss3)
confint(model_sss3)


#vigorous_cat_week 
model_sssss1 <- svyglm(Phenotypic_Age ~ vigorous_cat_week, data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssss1)
confint(model_sssss1)
model_sssss2 <- svyglm(Phenotypic_Age ~ vigorous_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssss2)
confint(model_sssss2)
model_sssss3 <- svyglm(Phenotypic_Age ~ vigorous_cat_week + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssss3)
confint(model_sssss3)

#vigorous_cat_day
model_sssssss1 <- svyglm(Phenotypic_Age ~ vigorous_cat_day, data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssssss1)
confint(model_sssssss1)
model_sssssss2 <- svyglm(Phenotypic_Age ~ vigorous_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssssss2)
confint(model_sssssss2)
model_sssssss3 <- svyglm(Phenotypic_Age ~ vigorous_cat_day + Gender + Race + CDAI + Marital_Status + Ratio_income_poverty
                  + BMI + sleep_disorders + Smoked_days + now_smoke + quit_smoking + Avg_alcohol_drinks + had_cancer
                  , data = Fulldat_Pheno, design = des, family = "gaussian") 
summary(model_sssssss3)
confint(model_sssssss3)