Last updated: 2022-10-23

Checks: 7 0

Knit directory: emlr_obs_preprocessing/

This reproducible R Markdown analysis was created with workflowr (version 1.7.0). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.


Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.

Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.

The command set.seed(20200707) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.

Great job! Recording the operating system, R version, and package versions is critical for reproducibility.

Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.

Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.

Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.

The results in this page were generated with repository version 221ee12. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.

Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:


Ignored files:
    Ignored:    .Rhistory
    Ignored:    .Rproj.user/
    Ignored:    data/
    Ignored:    output/

Untracked files:
    Untracked:  code/read_GLODAPv2_2022.Rmd

Unstaged changes:
    Modified:   code/Workflowr_project_managment.R

Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.


These are the previous versions of the repository in which changes were made to the R Markdown (analysis/read_Sabine_2004.Rmd) and HTML (docs/read_Sabine_2004.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.

File Version Author Date Message
html bafeecc jens-daniel-mueller 2022-06-07 Build site.
Rmd 46f2c6b jens-daniel-mueller 2022-06-06 included Nicos xover analysis for adjusted Knorr data
html b6bf005 jens-daniel-mueller 2022-04-26 Build site.
Rmd 5dfb9b0 jens-daniel-mueller 2022-04-26 checked unmasked column inventories of S04
html e949567 jens-daniel-mueller 2022-04-13 Build site.
html aea9afe jens-daniel-mueller 2022-04-07 Build site.
Rmd af08e38 jens-daniel-mueller 2022-04-07 rerun all with lat max 65N and without arcic
html f088f55 jens-daniel-mueller 2022-04-01 Build site.
html dde77eb jens-daniel-mueller 2022-04-01 Build site.
Rmd a1ea47d jens-daniel-mueller 2022-04-01 rerun all including arctic and North Atlantic biome
html 25c61af jens-daniel-mueller 2022-02-17 Build site.
Rmd 1c01415 jens-daniel-mueller 2022-02-17 added S04 budget calculation
html 6e65117 jens-daniel-mueller 2022-02-16 Build site.
html f2871b9 jens-daniel-mueller 2021-11-20 Build site.
html 0908ee5 jens-daniel-mueller 2021-11-15 Build site.
html e49875a jens-daniel-mueller 2021-07-07 Build site.
html 6312bd4 jens-daniel-mueller 2021-07-07 Build site.
Rmd 4905409 jens-daniel-mueller 2021-07-07 rerun with new setup_obs.Rmd file
html 58bc706 jens-daniel-mueller 2021-07-06 Build site.
Rmd 0db89e1 jens-daniel-mueller 2021-07-06 rerun with revised variable names

path_sabine_2004    <- "/nfs/kryo/work/updata/glodapv1_1/GLODAP_gridded.data/"
path_preprocessing  <- paste(path_root, "/observations/preprocessing/", sep = "")
library(marelac)

1 Data source

2 Read nc files

# read text files
AnthCO2_data <-
  read_csv(
    paste(path_sabine_2004,
          "AnthCO2.data/AnthCO2.data.txt",
          sep = ""),
    col_names = FALSE,
    na = "-999",
    col_types = list(.default = "d")
  )

# read respective depth layers and convert to vector
Depth_centers <-
  read_file(paste(path_sabine_2004,
                  "Depth.centers.txt",
                  sep = ""))

Depth_centers <- Depth_centers %>%
  str_split(",") %>%
  as_vector()

# read respective latitudes and convert to vector
Lat_centers <-
  read_file(paste(path_sabine_2004, "Lat.centers.txt",
                  sep = ""))

Lat_centers <- Lat_centers %>%
  str_split(",") %>%
  as_vector()

# read respective longitudes and convert to vector
Long_centers <-
  read_file(paste(path_sabine_2004, "Long.centers.txt",
                  sep = ""))

Long_centers <- Long_centers %>%
  str_split(",") %>%
  as_vector()

# match lon, lat and depth vectors with Cant value file
names(AnthCO2_data) <- Lat_centers

Long_Depth <-
  expand_grid(depth = Depth_centers, lon = Long_centers) %>%
  mutate(lon = as.numeric(lon),
         depth = as.numeric(depth))

tcant_3d <- bind_cols(AnthCO2_data, Long_Depth)

# adjust file dimensions
tcant_3d <- tcant_3d %>%
  pivot_longer(1:180, names_to = "lat", values_to = "tcant") %>%
  mutate(lat = as.numeric(lat))

tcant_3d <- tcant_3d %>%
  drop_na()

# harmonize coordinates
tcant_3d <- tcant_3d %>%
  mutate(lon = if_else(lon < 20, lon + 360, lon))

rm(AnthCO2_data,
   Long_Depth,
   Depth_centers,
   Lat_centers,
   Long_centers)

3 Apply basin mask

# use only three basin to assign general basin mask
# ie this is not specific to the MLR fitting
basinmask <- basinmask %>% 
  filter(MLR_basins == "2") %>% 
  select(lat, lon, basin_AIP)

tcant_3d_unmasked <- tcant_3d
tcant_3d <- inner_join(tcant_3d, basinmask)


ggplot() +
  geom_tile(data = tcant_3d_unmasked %>% 
              distinct(lon, lat),
            aes(lon, lat, fill = "basin mask not applied")) +
  geom_tile(data = tcant_3d %>% 
              distinct(lon, lat),
            aes(lon, lat, fill = "basin mask applied")) +
  coord_quickmap()

Version Author Date
b6bf005 jens-daniel-mueller 2022-04-26
rm(tcant_3d_unmasked)

4 Calculation

4.1 Column inventory

tcant_3d <- tcant_3d %>% 
  mutate(tcant_pos = if_else(tcant <= 0, 0, tcant))

tcant_inv_layers <- m_tcant_inv(tcant_3d)

tcant_inv <- tcant_inv_layers %>% 
  filter(inv_depth == params_global$inventory_depth_standard)

4.2 Zonal mean section

tcant_zonal <- m_zonal_mean_sd(tcant_3d)
m_dcant_budget(
  tcant_inv_layers %>%
    rename(dcant = tcant,
           dcant_pos = tcant_pos) %>%
    mutate(method = "total",
           data_source = "obs")) %>%
  select(-c(data_source, method)) %>% 
  group_by(estimate) %>%
  mutate(ratio = round(value / lag(value),3)) %>%
  ungroup() %>% 
  arrange(estimate, inv_depth)
# A tibble: 10 × 4
   inv_depth estimate  value  ratio
       <dbl> <chr>     <dbl>  <dbl>
 1       100 dcant      16.8 NA    
 2       500 dcant      63.0  3.75 
 3      1000 dcant      87.9  1.39 
 4      3000 dcant     102.   1.16 
 5     10000 dcant      97.7  0.963
 6       100 dcant_pos  16.8 NA    
 7       500 dcant_pos  63.0  3.75 
 8      1000 dcant_pos  88.2  1.40 
 9      3000 dcant_pos 104.   1.18 
10     10000 dcant_pos 106.   1.02 

5 Plots

5.1 Inventory map

p_map_cant_inv(
  df = tcant_inv,
  var = "tcant_pos",
  breaks = seq(0,max(tcant_inv$tcant_pos),5))

Version Author Date
bafeecc jens-daniel-mueller 2022-06-07
aea9afe jens-daniel-mueller 2022-04-07
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
0908ee5 jens-daniel-mueller 2021-11-15
58bc706 jens-daniel-mueller 2021-07-06
p_map_cant_inv(
  df = tcant_inv,
  var = "tcant",
  breaks = seq(0,max(tcant_inv$tcant_pos),5))

Version Author Date
bafeecc jens-daniel-mueller 2022-06-07
b6bf005 jens-daniel-mueller 2022-04-26

5.2 Horizontal plane maps

p_map_climatology(
  df = tcant_3d,
  var = "tcant",
  col = "divergent")

Version Author Date
aea9afe jens-daniel-mueller 2022-04-07
dde77eb jens-daniel-mueller 2022-04-01
0908ee5 jens-daniel-mueller 2021-11-15
58bc706 jens-daniel-mueller 2021-07-06

5.3 Global section

p_section_global(
  df = tcant_3d,
  var = "tcant",
  col = "divergent")

Version Author Date
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
58bc706 jens-daniel-mueller 2021-07-06

5.4 Sections at regular longitudes

p_section_climatology_regular(
  df = tcant_3d,
  var = "tcant",
  col = "divergent")

Version Author Date
aea9afe jens-daniel-mueller 2022-04-07
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
58bc706 jens-daniel-mueller 2021-07-06

5.5 Write files

tcant_3d %>%
  write_csv(paste(path_preprocessing,
                  "S04_tcant_3d.csv", sep = ""))

tcant_inv %>%
  write_csv(paste(path_preprocessing,
                  "S04_tcant_inv.csv", sep = ""))

tcant_zonal %>%
  write_csv(paste(path_preprocessing,
                  "S04_tcant_zonal.csv", sep = ""))

6 Anomalous changes

tcant_inv_S04 <- tcant_inv

dcant_inv_G19 <- read_csv(paste(path_preprocessing,
                               "G19_dcant_inv.csv", sep = ""))

7 Comparison of previous estimates

Cant inventory estimates of S04 (Sabine et al, 2004) and G19 (Gruber et al, 2019) were compared.

7.1 Merge data sets

cant_inv <- full_join(dcant_inv_G19 %>%
                        mutate(estimate = "G19") %>% 
                        rename(cant_pos = dcant_pos) %>% 
                        select(-dcant),
                      tcant_inv_S04 %>% 
                        mutate(estimate = "S04") %>% 
                        rename(cant_pos = tcant_pos) %>% 
                        select(-tcant))

rm(dcant_inv_G19, tcant_inv_S04)

7.2 Inventory maps

Spanning different time periods, the Cant inventories differ in magnitude. Please note, that we refer to cant_pos here, but strictly speaking we compare dcant and tcant.

map +
  geom_raster(data = cant_inv,
              aes(lon, lat, fill = cant_pos)) +
  scale_fill_viridis_c() +
  facet_wrap( ~ estimate, ncol = 1) +
  theme(
    axis.title = element_blank(),
    axis.text = element_blank(),
    axis.ticks = element_blank()
  )

Version Author Date
aea9afe jens-daniel-mueller 2022-04-07
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
6312bd4 jens-daniel-mueller 2021-07-07

7.3 Cant budgets

Global Cant inventories were estimated in Pg-C. Please note that here we only added positive Cant values in the upper m and do not apply additional corrections for areas not covered.

cant_inv <- cant_inv %>% 
  mutate(surface_area = earth_surf(lat, lon),
         cant_pos_grid = cant_pos*surface_area)

cant_inv_budget <- cant_inv %>% 
  group_by(estimate, basin_AIP) %>% 
  summarise(cant_pos_total = sum(cant_pos_grid)*12*1e-15,
            cant_pos_total = round(cant_pos_total,1)) %>% 
  ungroup() %>% 
  pivot_wider(values_from = cant_pos_total, names_from = basin_AIP) %>% 
  mutate(total = Atlantic + Indian + Pacific)

cant_inv_budget
# A tibble: 2 × 5
  estimate Atlantic Indian Pacific total
  <chr>       <dbl>  <dbl>   <dbl> <dbl>
1 G19          11      7.1    13.4  31.5
2 S04          39.6   23.4    41.4 104. 

7.4 Relative inventories

cant_inv_wide <- cant_inv %>%
  pivot_wider(values_from = c(cant_pos, cant_pos_grid),
              names_from = estimate)

cant_inv_wide <- cant_inv_wide %>% 
  drop_na() %>% 
  mutate(G19_rel = cant_pos_grid_G19 / sum(cant_pos_grid_G19),
         S04_rel = cant_pos_grid_S04 / sum(cant_pos_grid_S04),
         cant_ratio_rel = G19_rel / S04_rel)

cant_inv_rel <- cant_inv_wide %>%
  pivot_longer(
    cols = c(G19_rel, S04_rel),
    names_to = "estimate",
    values_to = "cant_pos_rel"
  )
map +
  geom_raster(data = cant_inv_rel,
              aes(lon, lat, fill = cant_pos_rel*100)) +
  scale_fill_viridis_c() +
  facet_wrap( ~ estimate, ncol = 1) +
  theme(
    axis.title = element_blank(),
    axis.text = element_blank(),
    axis.ticks = element_blank()
  )

Version Author Date
aea9afe jens-daniel-mueller 2022-04-07
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
6312bd4 jens-daniel-mueller 2021-07-07

7.5 Relative inventory ratios

map +
  geom_contour_filled(data = cant_inv_wide %>%
  filter(cant_ratio_rel < 10,
         cant_ratio_rel > 0.1),
                      aes(lon, lat, z = log10(cant_ratio_rel))) +
  coord_quickmap(expand = 0) +
  scale_fill_brewer(palette = "RdBu", direction = -1) +
  labs(title = "Cant inventory distribution | 1994-2007 vs preind-1994",
       subtitle = "Log ratio of relative contributions to total inventory") +
  theme(
    axis.title = element_blank(),
    axis.text = element_blank(),
    axis.ticks = element_blank(),
    legend.title = element_blank()
  )

Version Author Date
aea9afe jens-daniel-mueller 2022-04-07
f088f55 jens-daniel-mueller 2022-04-01
dde77eb jens-daniel-mueller 2022-04-01
6312bd4 jens-daniel-mueller 2021-07-07

sessionInfo()
R version 4.1.2 (2021-11-01)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: openSUSE Leap 15.3

Matrix products: default
BLAS:   /usr/local/R-4.1.2/lib64/R/lib/libRblas.so
LAPACK: /usr/local/R-4.1.2/lib64/R/lib/libRlapack.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] geomtextpath_0.1.0 colorspace_2.0-2   marelac_2.1.10     shape_1.4.6       
 [5] ggforce_0.3.3      metR_0.11.0        scico_1.3.0        patchwork_1.1.1   
 [9] collapse_1.7.0     forcats_0.5.1      stringr_1.4.0      dplyr_1.0.7       
[13] purrr_0.3.4        readr_2.1.1        tidyr_1.1.4        tibble_3.1.6      
[17] ggplot2_3.3.5      tidyverse_1.3.1    workflowr_1.7.0   

loaded via a namespace (and not attached):
 [1] fs_1.5.2           bit64_4.0.5        lubridate_1.8.0    gsw_1.0-6         
 [5] RColorBrewer_1.1-2 httr_1.4.2         rprojroot_2.0.2    tools_4.1.2       
 [9] backports_1.4.1    bslib_0.3.1        utf8_1.2.2         R6_2.5.1          
[13] DBI_1.1.2          withr_2.4.3        tidyselect_1.1.1   processx_3.5.2    
[17] bit_4.0.4          compiler_4.1.2     git2r_0.29.0       textshaping_0.3.6 
[21] cli_3.1.1          rvest_1.0.2        xml2_1.3.3         isoband_0.2.5     
[25] labeling_0.4.2     sass_0.4.0         scales_1.1.1       checkmate_2.0.0   
[29] SolveSAPHE_2.1.0   callr_3.7.0        systemfonts_1.0.3  digest_0.6.29     
[33] rmarkdown_2.11     oce_1.5-0          pkgconfig_2.0.3    htmltools_0.5.2   
[37] highr_0.9          dbplyr_2.1.1       fastmap_1.1.0      rlang_1.0.2       
[41] readxl_1.3.1       rstudioapi_0.13    jquerylib_0.1.4    generics_0.1.1    
[45] farver_2.1.0       jsonlite_1.7.3     vroom_1.5.7        magrittr_2.0.1    
[49] Rcpp_1.0.8         munsell_0.5.0      fansi_1.0.2        lifecycle_1.0.1   
[53] stringi_1.7.6      whisker_0.4        yaml_2.2.1         MASS_7.3-55       
[57] grid_4.1.2         parallel_4.1.2     promises_1.2.0.1   crayon_1.4.2      
[61] haven_2.4.3        hms_1.1.1          seacarb_3.3.0      knitr_1.37        
[65] ps_1.6.0           pillar_1.6.4       reprex_2.0.1       glue_1.6.0        
[69] evaluate_0.14      getPass_0.2-2      data.table_1.14.2  modelr_0.1.8      
[73] vctrs_0.3.8        tzdb_0.2.0         tweenr_1.0.2       httpuv_1.6.5      
[77] cellranger_1.1.0   gtable_0.3.0       polyclip_1.10-0    assertthat_0.2.1  
[81] xfun_0.29          broom_0.7.11       later_1.3.0        viridisLite_0.4.0 
[85] ellipsis_0.3.2