#### 1_set_up_n_create_time_activity_diagram_individual_n_aggregate.R #### library(tidyverse) library(sp) library(maptools) library(templates) library(raster) library(rgdal) library(haven) library(data.table) library(car) library(RColorBrewer) library(varhandle) library(grDevices) library(dplyr) library(plotrix) library(graphics) #### PART ONE: Set up parameters, load and process data #### # Set up work directory setwd("C:\\Timeline_diagrams") getwd() ## set parameters # deprivation group myDep <- 1 # width factor of % activity actAgglwd <- 0.9 # width factor of individual timelines actIndlwd <- 2 myDiffFactor = 9 myrds_name <- paste("IMD", myDep, "_obese_", "active_Boy_",sep='') myLefend2 <- 110 ## Set X, Y #### Xleft <- -600 Xright <- 3600 + Xleft # center of aggregates for left and right panals Xcl = -300 Xcr = Xright - 300 # set Y at bottom and top Yb <- 6*60*2 Yt <- 24*60*2 ## set color #### # dark green for MVPA col4 <- "#33a02c" col4a <- "#1b7837CC" # light green for LPA col2 <- "#a6dba0" col2a <- "#a6dba0CC" # orange for SB col1 <- "#fdbf6f" col1a <- "#fdbf6fBF" # purple for sleep col5 <- "#c2a5cf" # set grey and black color col_grey <- "#525252" col_lightgrey <- "#d9d9d9" col_black <- "#252525" col_nodata <- "#f0f0f0" ## Load datasets # read hour (06:00 to 24:00) dataset df.hour <- readRDS("df.hour.rds") # read minute dataset (at 10 minute interval) df.t_num <- readRDS("df.t_num.rds") # read preprocessed accelerometer data df_all.uniq_6_24_depq1_boy <- readRDS("df_all.uniq_6_24_depq1_boy.rds") # Filter dataset to Wednesday and Thursday and distance between home and school <= 3000m df_all.uniq_6_24_depq1_boy <- df_all.uniq_6_24_depq1_boy %>% filter((DAY == "Wednesday" | DAY == "Thursday") & df_all.uniq_6_24_depq1_boy$H_Shc_Dist < 3000) # identify valid days with >= 7 hour PAs and Sed time #### # summarise time and activity sumpids <- df_all.uniq_6_24_depq1_boy %>% filter(actnew != 0, actnew != 5) %>% group_by(pid, date_num, episode) %>% summarise(sum_pid_date= n()) %>% ungroup() %>% group_by(pid, date_num) %>% summarise(sum_acts= sum(sum_pid_date), sum_episode = n()) %>% ungroup() # Filter data to those with more than 7 hours (7*60*2= 840) spent on PA and SB sumpids2 <- sumpids %>% filter((sum_acts - sum_episode) >=840) %>% mutate(pid_date = paste0(pid,date_num)) unique(sumpids2$pid_date) # Get unique pid sumpids3 <- sumpids2 %>% dplyr::select(pid_date) %>% mutate(id = row_number()) # Get start position of individual timelines no.boys <- length(unique(sumpids3$id)) xposi <- (Xright - Xleft - 100 - 20 - 1200) / no.boys # filter data to valid days df_all.uniq_6_24_depq1_boy <- df_all.uniq_6_24_depq1_boy %>% mutate(pid_date = paste0(pid,date_num)) df_valid_6_24_dep1_boy <- df_all.uniq_6_24_depq1_boy %>% inner_join(sumpids3, by="pid_date") %>% distinct() # change ovwt_obese to two groups df_valid_6_24_dep1_boy <- df_valid_6_24_dep1_boy %>% mutate(ovwt_obese = ifelse(ovwt_obese <= 2, 1, 3)) unique(df_valid_6_24_dep1_boy$ovwt_obese) saveRDS(df_valid_6_24_dep1_boy, "df_active_6_24_dep1_boy.rds") # Reduce memoery rm(sumpids) rm(sumpids2) rm(sumpids3) rm(df_all.uniq_6_24_depq1_boy) gc() memory.size() # get order ID by obesity status, travel mode, home toschool distance, and id. boy_active_plus0 <- df_valid_6_24_dep1_boy %>% dplyr::select("ovwt_obese","STActive", "H_Shc_Dist", "id") boy_active_plus_order <- boy_active_plus0 %>% distinct(id, STActive, H_Shc_Dist, ovwt_obese, keep_all = FALSE) # sort boy_active_plus_order <- boy_active_plus_order %>% arrange(ovwt_obese, STActive, H_Shc_Dist,id) # get X poositions for individual timelines boy_active_plus_orderID <- boy_active_plus_order %>% tibble::rowid_to_column("X_Obese_Dist") unique(boy_active_plus_orderID$ovwt_obese) unique(boy_active_plus_orderID$STActive) unique(boy_active_plus_orderID$X_Obese_Dist) # link order id to data and get X field boy_active_plus_XObeseDist <- df_valid_6_24_dep1_boy %>% left_join(boy_active_plus_orderID %>% dplyr::select(-STActive, -ovwt_obese, -H_Shc_Dist), by="id") %>% distinct() %>% mutate(X_Obese_ActiveT = X_Obese_Dist*xposi + (xposi / 2)) unique(boy_active_plus_XObeseDist$X_Obese_ActiveT) # Reduce memoery rm(boy_active_plus0) rm(boy_active_plus_order) rm(boy_active_plus_orderID) gc() memory.size() # add space between normal and obese groups boy_active_plus_XObeseDist <- boy_active_plus_XObeseDist %>% mutate(X_Obese_Dist = ifelse(ovwt_obese <= 2, X_Obese_ActiveT, X_Obese_ActiveT+100), col = ifelse(actnew == 0, "#f0f0f0", ifelse(actnew == 1, "#fdbf6f", ifelse(actnew == 2, "#a6dba0", ifelse(actnew == 4, "#1b7837", ifelse(actnew == 5, "#c2a5cf", 0)))))) # link data to minute of time boy_active_plus_XObeseDist_t_num <- boy_active_plus_XObeseDist %>% left_join(df.t_num, by="time_num") %>% distinct() %>% drop_na(min2) n.boy <- length(boy_active_plus_XObeseDist_t_num$id) # mergeBoy# boy_active_plus_XObeseDist_path <- boy_active_plus_XObeseDist %>% distinct(pid, id, ovwt_obese, actnew, episode, col, X_Obese_Dist, keep_all = TRUE) ## add color #### boy_active_plus_XObeseDist_path$col2 <- '' boy_active_plus_XObeseDist_path$col2 <- ifelse(boy_active_plus_XObeseDist_path$actnew == 0, "#f0f0f0", ifelse(boy_active_plus_XObeseDist_path$actnew == 1, "#fdbf6f", ifelse(boy_active_plus_XObeseDist_path$actnew == 2, "#a6dba0", ifelse(boy_active_plus_XObeseDist_path$actnew == 4, "#1b7837", ifelse(boy_active_plus_XObeseDist_path$actnew == 5, "#c2a5cf", 0))))) # create visualisation using spatial approach # # create spatial points coordinates(boy_active_plus_XObeseDist)=~ X_Obese_Dist+time_num # this part is slow x <- lapply(split(boy_active_plus_XObeseDist, boy_active_plus_XObeseDist$episode), function(x) Lines(list(Line(coordinates(x))), x$episode[1L])) # use parallel programming to speed up a little bit. # temp <- split(boy_active_plus_XObeseDist, # boy_active_plus_XObeseDist$episode) # # library(parallel) # cl <- makeCluster(12) # # clusterEvalQ(cl, require(sp)) # x <- parLapply(cl, temp, function(x) # Lines(Line(coordinates(x)), # x$episode[1L])) # # stopCluster(cl) lines_boy_active_plus_XObeseDist <- SpatialLines(x) #sort boy_active_plus_XObeseDist_path <- boy_active_plus_XObeseDist_path %>% arrange(pid, episode) %>% as.data.frame() # add data to SpatialLinesDataFrame rownames(boy_active_plus_XObeseDist_path) <- boy_active_plus_XObeseDist_path$episode spdflines_boy_active_plus_XObeseDist <- SpatialLinesDataFrame(lines_boy_active_plus_XObeseDist, data.frame(data =boy_active_plus_XObeseDist_path), match.ID=TRUE) rm(lines_boy_active_plus_XObeseDist) #### normal Boys #### # calculate no. of normal Boys mergeBoy0 <- boy_active_plus_XObeseDist_t_num %>% filter(ovwt_obese <= 2) n.Boy0 <- length(unique(mergeBoy0$id)) n.Boy0 # get total count of act by min and act group mergeBoygrp <- mergeBoy0 %>% filter(actnew != 0 & !is.na(min2)) %>% group_by(min2, actnew) %>% summarise(n_act= n()) sumacts <- mergeBoygrp %>% group_by(min2) %>% summarise(sumacts= sum(n_act)) sumnormal <- mergeBoygrp %>% left_join(sumacts, by="min2") %>% distinct() %>% mutate(mean_act_pct = n_act/sumacts*100) BoyDistMin <- mergeBoy0 %>% pull(H_Shc_Dist) %>% min %>% round(0) %>% as.character() BoyDistMax <- mergeBoy0 %>% pull(H_Shc_Dist) %>% max %>% round(0) %>% as.character() # create aggregate for act 5#### sumnormal5 <- sumnormal %>% filter(actnew==5) sumnormal5X <- sumnormal5 %>% left_join(df.t_num, by="min2") %>% distinct() # sumnormal5X$min2 <- NULL sumnormal5X_1 <- sumnormal5X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumnormal5X2 <- sumnormal5X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col = ifelse(!is.na(mean_act_pct), "#c2a5cf", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA) )) sumnormal5X$X5 <- Xcl #### write sumnormal5x to rds #### saveRDS(sumnormal5X, paste(myrds_name,"sumnormal5X.rds",sep="")) saveRDS(sumnormal5X2, paste(myrds_name,"sumnormal5X2.rds",sep="")) coordinates(sumnormal5X)= ~X5 + time_num x <- lapply(split(sumnormal5X, sumnormal5X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumnormal5X <- SpatialLines(x) sumnormal5X2 <- as.data.frame(sumnormal5X2) # add data to SpatialLinesDataFrame rownames(sumnormal5X2) <- sumnormal5X2$min2 spdflines_sumnormal5X <- SpatialLinesDataFrame(lines_sumnormal5X, data.frame(data = sumnormal5X2), match.ID=TRUE) # create aggregate for act 1#### sumnormal1 <- sumnormal %>% filter(actnew==1) sumnormal1X <- sumnormal1 %>% left_join(df.t_num, by="min2") %>% distinct() sumnormal1X_1 <- sumnormal1X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumnormal1X2 <- sumnormal1X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#fdbf6f80", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumnormal1X$X1 <- Xcl #### write sumnormal1x to rds #### saveRDS(sumnormal1X, paste(myrds_name,"sumnormal1X.rds",sep="")) saveRDS(sumnormal1X2, paste(myrds_name,"sumnormal1X2.rds",sep="")) coordinates(sumnormal1X)= ~X1 + time_num x <- lapply(split(sumnormal1X, sumnormal1X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumnormal1X <- SpatialLines(x) sumnormal1X2 <- as.data.frame(sumnormal1X2) # add data to SpatialLinesDataFrame rownames(sumnormal1X2) <- sumnormal1X2$min2 spdflines_sumnormal1X <- SpatialLinesDataFrame(lines_sumnormal1X, data.frame(data =sumnormal1X2), match.ID=TRUE) # create aggregate for act 2#### sumnormal2 <- sumnormal %>% filter(actnew==2) sumnormal2X <- sumnormal2 %>% left_join(df.t_num, by="min2") %>% distinct() # sumnormal2X$min2 <- NULL sumnormal2X_1 <- sumnormal2X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumnormal2X2 <- sumnormal2X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#a6dba080", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumnormal2X$X2 <- Xcl #### write sumnormal2X to rds #### saveRDS(sumnormal2X, paste(myrds_name,"sumnormal2X.rds",sep="")) saveRDS(sumnormal2X2, paste(myrds_name,"sumnormal2X2.rds",sep="")) coordinates(sumnormal2X)= ~X2 + time_num x <- lapply(split(sumnormal2X, sumnormal2X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumnormal2X <- SpatialLines(x) sumnormal2X2 <- as.data.frame(sumnormal2X2) # add data to SpatialLinesDataFrame rownames(sumnormal2X2) <- sumnormal2X2$min2 spdflines_sumnormal2X <- SpatialLinesDataFrame(lines_sumnormal2X, data.frame(data =sumnormal2X2), match.ID=TRUE) # create aggregate for act 4#### sumnormal4 <- sumnormal %>% filter(actnew==4) sumnormal4X <- sumnormal4 %>% left_join(df.t_num, by="min2") %>% distinct() sumnormal4X_1 <- sumnormal4X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumnormal4X2 <- sumnormal4X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#1b7837BF", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumnormal4X$X4 <- Xcl #### write sumnormal4X to rds #### saveRDS(sumnormal4X, paste(myrds_name,"sumnormal4X.rds",sep="")) saveRDS(sumnormal4X2, paste(myrds_name,"sumnormal4X2.rds",sep="")) coordinates(sumnormal4X)= ~X4 + time_num x <- lapply(split(sumnormal4X, sumnormal4X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumnormal4X <- SpatialLines(x) sumnormal4X2 <- as.data.frame(sumnormal4X2) # add data to SpatialLinesDataFrame rownames(sumnormal4X2) <- sumnormal4X2$min2 spdflines_sumnormal4X <- SpatialLinesDataFrame(lines_sumnormal4X, data.frame(data =sumnormal4X2), match.ID=TRUE) #### Boy with obese #### # calculate no of obese Boys mergeBoy1 <- boy_active_plus_XObeseDist_t_num %>% filter(ovwt_obese >= 3) n.Boy1 <- length(unique(mergeBoy1$id)) n.Boy1 # get total count of act by min and act group mergeBoygrp <- mergeBoy1 %>% filter(actnew != 0 & !is.na(min2))%>% group_by(min2, actnew) %>% summarise(n_act= n()) sumacts <- mergeBoygrp %>% group_by(min2) %>% summarise(sumacts= sum(n_act)) sumobese <- mergeBoygrp %>% left_join(sumacts, by="min2") %>% distinct() %>% mutate(mean_act_pct = (n_act/sumacts)*100) BoyDistMinobese <- mergeBoy1 %>% pull(H_Shc_Dist) %>% min %>% round(0) %>% as.character() BoyDistMaxobese <- mergeBoy1 %>% pull(H_Shc_Dist) %>% max %>% round(0) %>% as.character() # create aggregate for act 5#### sumobese5 <- sumobese %>% filter(actnew==5) sumobese5X <- sumobese5 %>% left_join(df.t_num, by="min2") %>% distinct() sumobese5X_1 <- sumobese5X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumobese5X2 <- sumobese5X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#c2a5cf", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumobese5X$X5 <- Xcr #### write sumobese5x to rds #### saveRDS(sumobese5X, paste(myrds_name,"sumobese5X.rds",sep="")) saveRDS(sumobese5X2, paste(myrds_name,"sumobese5X2.rds",sep="")) coordinates(sumobese5X)= ~X5 + time_num x <- lapply(split(sumobese5X, sumobese5X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumobese5X <- SpatialLines(x) sumobese5X2 <- as.data.frame(sumobese5X2) # add data to SpatialLinesDataFrame rownames(sumobese5X2) <- sumobese5X2$min2 spdflines_sumobese5X <- SpatialLinesDataFrame(lines_sumobese5X, data.frame(data =sumobese5X2), match.ID=TRUE) # create aggregate for act 1#### sumobese1 <- sumobese %>% filter(actnew==1) sumobese1X <- sumobese1 %>% left_join(df.t_num, by="min2") %>% distinct() sumobese1X_1 <- sumobese1X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumobese1X2 <- sumobese1X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#fdbf6f80", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumobese1X$X1 <- Xcr #### write sumobese1x to rds #### saveRDS(sumobese1X, paste(myrds_name,"sumobese1X.rds",sep="")) saveRDS(sumobese1X2, paste(myrds_name,"sumobese1X2.rds",sep="")) coordinates(sumobese1X)= ~X1 + time_num x <- lapply(split(sumobese1X, sumobese1X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumobese1X <- SpatialLines(x) sumobese1X2 <- as.data.frame(sumobese1X2) # add data to SpatialLinesDataFrame rownames(sumobese1X2) <- sumobese1X2$min2 spdflines_sumobese1X <- SpatialLinesDataFrame(lines_sumobese1X, data.frame(data =sumobese1X2), match.ID=TRUE) # create aggregate for act 2#### sumobese2<- sumobese %>% filter(actnew==2) sumobese2X <- sumobese2 %>% left_join(df.t_num, by="min2") %>% distinct() sumobese2X_1 <- sumobese2X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumobese2X2 <- sumobese2X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#a6dba080", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumobese2X$X2 <- Xcr #### write sumobese2X to rds #### saveRDS(sumobese2X, paste(myrds_name,"sumobese2X.rds",sep="")) saveRDS(sumobese2X2, paste(myrds_name,"sumobese2X2.rds",sep="")) coordinates(sumobese2X)= ~X2 + time_num x <- lapply(split(sumobese2X, sumobese2X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumobese2X <- SpatialLines(x) sumobese2X2 <- as.data.frame(sumobese2X2) # add data to SpatialLinesDataFrame rownames(sumobese2X2) <- sumobese2X2$min2 spdflines_sumobese2X <- SpatialLinesDataFrame(lines_sumobese2X, data.frame(data =sumobese2X2), match.ID=TRUE) # create aggregate for act 4#### sumobese4 <- sumobese %>% filter(actnew==4) sumobese4X <- sumobese4 %>% left_join(df.t_num, by="min2") %>% distinct() sumobese4X_1 <- sumobese4X %>% distinct(actnew, n_act, sumacts, min2, time_num, keep_all = FALSE) sumobese4X2 <- sumobese4X %>% distinct(mean_act_pct, min2, keep_all = FALSE) %>% mutate(col =ifelse(!is.na(mean_act_pct), "#1b7837BF", ifelse((is.na(mean_act_pct) | mean_act_pct == 0), NA, NA))) sumobese4X$X4 <- Xcr #### write sumobese4X to rds #### saveRDS(sumobese4X, paste(myrds_name,"sumobese4X.rds",sep="")) saveRDS(sumobese4X2, paste(myrds_name,"sumobese4X2.rds",sep="")) coordinates(sumobese4X)= ~X4 + time_num x <- lapply(split(sumobese4X, sumobese4X$min2), function(x) Lines(list(Line(coordinates(x))), x$min2[1L])) lines_sumobese4X <- SpatialLines(x) sumobese4X2 <- as.data.frame(sumobese4X2) # add data to SpatialLinesDataFrame rownames(sumobese4X2) <- sumobese4X2$min2 spdflines_sumobese4X <- SpatialLinesDataFrame(lines_sumobese4X, data.frame(data =sumobese4X2), match.ID=TRUE) rm(boy_active_plus_XObeseDist) rm(boy_active_plus_XObeseDist_t_num) gc() memory.size() #### PART TWO: Plot Time-Activity Diagram showing both individual and aggregated information #### # set file name for plot #### myfile<-paste("Obese_Boys_IMDQ", myDep, "_individual_n_aggregate_test2.png",sep='') png(myfile, width=21, height=15, units = "cm", res = 2400, pointsize = 10) # add an empty plot #### plot(1, type="n", xlab="", ylab="", xlim=c((Xleft - 100), (Xright+20)), ylim=c((Yb-320), (Yt + 120)), axes=FALSE, ann=FALSE) # draw rectangle rect((n.Boy0*xposi+105), 6*120, (Xright - 600 - 5), 24*120, density = NULL, col = col_nodata, border = NULL) rect(((Xright - 600 -5)), 6*120, Xright, 24*120, density = NULL, col = "white", border = NULL) rect(0, 6*120, (n.Boy0*xposi + 15), 24*120, density = NULL, col = col_nodata, border = NULL) rect(Xleft, 6*120, 5, 24*120, density = NULL, col = "white", border = NULL) ## plot hour referencing lines for (i in seq(6, 24, by = 3)) { # hour referencing lines segments((Xleft), i*120, Xright, i*120, col = col_lightgrey, lty = , lwd = 0.4) #label hour lab <- paste(c(i, ":00"), collapse = "") text((Xleft), i*120, lab, cex = .9, pos = 2) } ## plot individual activities #### plot(spdflines_boy_active_plus_XObeseDist, lwd =actIndlwd, col = boy_active_plus_XObeseDist_path$col, lend = 1, add = TRUE) rm(spdflines_boy_active_plus_XObeseDist) memory.size() gc() memory.size() ## plot aggregated activities for those classified as normal ## # plot aggregates 5 #### plot( spdflines_sumnormal5X, lwd = (data.frame(sumnormal5X2)$mean_act_pct*actAgglwd), col = col5, lend = 1, add = TRUE ) # plot aggregates 1 #### plot( spdflines_sumnormal1X, lwd = (data.frame(sumnormal1X2)$mean_act_pct*actAgglwd), col = (data.frame(sumnormal1X2)$col), lend = 1, add = TRUE ) # plot aggregates 2 #### plot( spdflines_sumnormal2X, lwd = (data.frame(sumnormal2X2)$mean_act_pct*actAgglwd), col = (data.frame(sumnormal2X2)$col), lend = 1, add = TRUE ) # plot aggregates 4 #### plot( spdflines_sumnormal4X, lwd = (data.frame(sumnormal4X2)$mean_act_pct*actAgglwd), col = (data.frame(sumnormal4X2)$col), lend = 1, add = TRUE ) ## plot aggregated activities for those classified as obese## # plot aggregates 5 #### plot( spdflines_sumobese5X, lwd = (data.frame(sumobese5X2)$mean_act_pct*actAgglwd), col = col5, lend = 1, add = TRUE ) # plot aggregates 1 #### plot( spdflines_sumobese1X, lwd = (data.frame(sumobese1X2)$mean_act_pct*actAgglwd), col = (data.frame(sumobese1X2)$col), lend = 1, add = TRUE ) # plot aggregates 2 #### plot( spdflines_sumobese2X, lwd = (data.frame(sumobese2X2)$mean_act_pct*actAgglwd), col = (data.frame(sumobese2X2)$col), lend = 1, add = TRUE ) # plot aggregates 4 #### plot( spdflines_sumobese4X, lwd = (data.frame(sumobese4X2)$mean_act_pct*actAgglwd), col = (data.frame(sumobese4X2)$col), lend = 1, add = TRUE ) #### legend #### pre_legend_title = "Weight status (classified using Cole 2012 IOTF cut-offs)" # add title on top legend("top", legend = "", title = pre_legend_title, cex=1.1, bty = "n", inset = -.01) # text((-30), (Yb-60), lab=as.factor(BoyDistMin), cex = .9, pos = 4, font=1) text((n.Boy0*xposi+60), (Yb-60), lab=BoyDistMax, cex = .9, pos = 2, font=1) text((n.Boy0*xposi + 60), (Yb-60), lab=BoyDistMinobese, cex = .9, pos = 4, font=1) text((n.Boy0*xposi + 150 + n.Boy1*xposi), (Yb-60), lab=BoyDistMaxobese, cex = .9, pos = 2, font=1) text((Xleft-30), (Yt+60), lab="Underweight or Normal", cex = 1.0, pos = 4, font=1) text((n.Boy0*xposi + 70), (Yt+60), lab="Overweight or Obese", cex = 1.0, pos = 4, font=1) text(((n.Boy0*xposi + 150 + n.Boy1*xposi) -(-30))/2, (Yb-70), lab="Individuals ordered by the shortest home to school distance (m)", cex = 0.9, pos = 1, font=1) # horizontal box lines segments((Xleft), Yt, (n.Boy0*xposi + 15), Yt, col = "black", lty = , lwd = 1) segments((Xleft), Yb, (n.Boy0*xposi + 15), Yb, col = "black", lty = , lwd = 1) segments((Xright), Yt, (n.Boy0*xposi+105), Yt, col = "black", lty = , lwd = 1) segments((Xright), Yb, (n.Boy0*xposi+105), Yb, col = "black", lty = , lwd = 1) # vertical box lines segments((Xleft), Yb, (Xleft), Yt, col = "black", lty = , lwd = 1) segments((n.Boy0*xposi + 15), Yb, (n.Boy0*xposi + 15), Yt, col = "black", lty = , lwd = 1) segments((Xright), Yb, (Xright), Yt, col = "black", lty = , lwd = 1) segments((n.Boy0*xposi+105), Yb, (n.Boy0*xposi+105), Yt, col = "black", lty = , lwd = 1) # legends and labels for act % # right side segments((Xcr), (Yb-61.5), (Xcr), (Yb-63.5), col = "black", lty = , lend = 1, lwd = (100*actAgglwd)) segments((Xcr - 140-120), (Yb - 52.5), (Xcr - 140-120), (Yb - 72.5), col = "black", lty = , lwd = 1) segments((Xcr + 140 +120), (Yb - 52.5), (Xcr + 140 +120), (Yb - 72.5), col = "black", lty = , lwd = 1) boxed.labels((Xcr), (Yb-62.5),"100%",bg="white",border=NA, cex = .9) text((Xcr), (Yb - 70), lab="Mean % of activities", cex = .9, pos = 1, font=1) # left side segments((Xcl), (Yb-61.5), (Xcl), (Yb-63.5), col = "black", lty = , lend = 1, lwd = (100*actAgglwd)) segments((Xcl - 140-120), (Yb - 52.5), (Xcl - 140-120), (Yb - 72.5), col = "black", lty = , lwd = 1) segments((Xcl + 140 +120), (Yb - 52.5), (Xcl + 140 +120), (Yb - 72.5), col = "black", lty = , lwd = 1) boxed.labels((Xcl), (Yb-62.5),"100%",bg="white",border=NA, cex = .9) text((Xcl), (Yb-70), lab="Mean % of activities", cex = .9, pos = 1, font=1) # create customised point legend segdist <- ((Xright - Xleft) + 200)/5 # dark green for MVPA points((Xleft+20), (Yb-270), bg = "#1b7837BF", cex = 2, cex.lab = 1.0, pch = 22, col = NA) points((Xleft+ 20+myLefend2), (Yb-270), bg = "#1b7837", cex = 2, cex.lab = 1.0,pch = 22, col = NA) text((Xleft+20 +100), (Yb-270), lab=" MVPA", cex = 1.0, pos = 4, font=1) # light green for LPA points((Xleft + segdist), (Yb-270), bg = "#a6dba080", cex = 2, cex.lab = 1.0, pch = 22, col = NA) points((Xleft + segdist + myLefend2), (Yb-270), bg = "#a6dba0", cex = 2, cex.lab = 1.0, pch = 22, col = NA) text((Xleft + segdist + 100), (Yb-270), lab=" LPA", cex = 1.0, pos = 4, font=1) # orange for SB points((Xleft + segdist*2), (Yb-270), bg = "#fdbf6f80", cex = 2, cex.lab = 1.0, pch = 22, col = NA) points((Xleft + segdist*2 + myLefend2), (Yb-270), bg = "#fdbf6f", cex = 2, cex.lab = 1.0, pch = 22, col = NA) text((Xleft + segdist*2 + 100), (Yb-270), lab=" SB", cex = 1.0, pos = 4, font=1) # purple for sleep points((Xleft + segdist*3), (Yb-270), bg = "#c2a5cf", cex = 2, cex.lab = 1.0, pch = 22, col = NA) points((Xleft + segdist*3 + myLefend2), (Yb-270), bg = "#c2a5cf", cex = 2, cex.lab = 1.0, pch = 22, col = NA) text((Xleft + segdist*3 + 100), (Yb-270), lab=" Sleep", cex = 1.0, pos = 4, font=1) # grey for no data points((Xleft + segdist*4), (Yb-270), bg = "#f0f0f0", cex = 2, cex.lab = 1.0, pch = 22, col = NA) text((Xleft + segdist*4 + 20), (Yb-270), lab=" No data", cex = 1.0, pos = 4, font=1) dev.off()