How can plot my own data in a grid in a map sf but return vacum

Question

I am trying to summarize some statistics in the grid that I made, however something fails when I try to do it.

This is my data

head(catk)
Simple feature collection with 6 features and 40 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 303.22 ymin: -61.43 xmax: 303.95 ymax: -60.78
Geodetic CRS:  WGS 84
# A tibble: 6 × 41
     X1   day month  year c1_id greenweight_caught_kg obs_haul_id
  <dbl> <dbl> <dbl> <dbl> <dbl>                 <dbl> <lgl>      
1     1     4    12  1997 26529                  7260 NA         
2     2     4    12  1997 26530                  7920 NA         
3     3     4    12  1997 26531                  4692 NA         
4     4     4    12  1997 26532                  5896 NA         
5     5     4    12  1997 26533                    88 NA         
6     6     5    12  1997 26534                  7040 NA         
# … with 34 more variables: obs_logbook_id <lgl>, obs_haul_number <lgl>,
#   haul_number <dbl>, vessel_name <chr>, vessel_nationality_code <chr>,
#   fishing_purpose_code <chr>, season_ccamlr <dbl>,
#   target_species <chr>, asd_code <dbl>, trawl_technique <lgl>,
#   catchperiod_code <chr>, date_catchperiod_start <date>,
#   datetime_set_start <dttm>, datetime_set_end <dttm>,
#   datetime_haul_start <dttm>, datetime_haul_end <dttm>, …

and I did this raster

an <- getData("GADM", country = "ATA", level = 0)
an@data$NAME_0

e <- extent(-70,-40,-68,-60)
rc <- crop(an, e)
proj4string(rc) <- CRS("+init=epsg:4326")
rc3 <- st_as_sf(rc)

catk <- st_as_sf(catk, coords = c("Longitude", "Latitude"),  crs = 4326) %>%
  st_shift_longitude()

Grid <- rc3 %>%
  st_make_grid(cellsize = c(1,0.4)) %>% # para que quede cuadrada
  st_cast("MULTIPOLYGON") %>%
  st_sf() %>%
  mutate(cellid = row_number()) 

result <- Grid %>%
  st_join(catk) %>% 
  group_by(cellid) %>% 
  summarise(sum_cat = sum(Catcht))

but I can not represent the data in the grid

ggplot() +
  geom_sf(data = Grid, color="#d9d9d9", fill=NA) +
  geom_sf(data = rc3) +
  theme_bw() + 
  coord_sf() +
  scale_alpha(guide="none")+
  xlab(expression(paste(Longitude^o,~'O'))) +
  ylab(expression(paste(Latitude^o,~'S')))+
  guides( colour = guide_legend()) +
  theme(panel.background = element_rect(fill = "#f7fbff"),
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank())+
  theme(legend.position = "right")+
  xlim(-69,-45)

fail plot

Please help me to find this solution!!

Hi @MauroMardones. Could you please share `catk` object or at least the initial file that you read it from? Another option is to run `dput(head(catk ))` and share the output ` — dieghernan, Mar 22 '22 at 07:52
Hi @dieghernan. I read a .csv file and it called `catk`. ``` class(catk) [1] "spec_tbl_df" "tbl_df" "tbl" "data.frame" ``` Thanks in advice! — MauroMardones, Mar 22 '22 at 14:58
Thanks @MauroMardones, what I meant is to get the object (or even better, the csv), not the type of object. It is very hard to figure out how to helkp you without access to your data — dieghernan, Mar 22 '22 at 15:40
I understand, but how can post a large *.csv* file?? (sorry) — MauroMardones, Mar 22 '22 at 18:06
Just a few lines (~20) are enough, if you have GitHub you can upload it into a gist and copy the link. Another option is, as I said, run `dput(head(catk ))` right after you create the object and update your question pasting the output of the console — dieghernan, Mar 22 '22 at 19:20
See an example https://stackoverflow.com/questions/66724987/the-output-of-dputheaddata-20-in-data-frame-in-r — dieghernan, Mar 22 '22 at 19:21
You can also check the `reprex` package about some guidelines on creating a bold **repr**oducible **ex**ample https://reprex.tidyverse.org/articles/reprex-dos-and-donts.html — dieghernan, Mar 22 '22 at 19:26

score 0 · Accepted Answer · answered Mar 23 '22 at 08:09

So I just saw that you shifted the coordinates with st_shift_longitude() and therefore your bounding box is:

Bounding box:  xmin: 303.22 ymin: -61.43 xmax: 303.95 ymax: -60.78

Do you really need it? That doesn't match with your defined extent

e <- extent(-70,-40,-68,-60)

And a bbox for WGS84 is suppose to be at max c(-180, -90, 180, 90).

Also, on your plot you are not instructed ggplot2 to do anything with the values of catk. Grid and rc3 do not have anything from catk, is the result object.

Anyway, I give a try to your problem even though I don't have access to your dataset. I represent on each cell sum_cat from result

library(raster)
library(sf)
library(dplyr)
library(ggplot2)

# Mock your data
catk <- structure(list(Longitude = c(-59.0860203764828, -50.1352159580643, 
                                     -53.7671292009259, -67.9105254106185, -67.5753491797446, -51.7045571975837, 
                                     -45.6203830411619, -61.2695183762776, -51.6287384188695, -52.244074640978, 
                                     -45.4625770258213, -51.0935832496694, -45.6375681312716, -44.744215508174, 
                                     -66.3625310507564), Latitude = c(-62.0038884948778, -65.307178606448, 
                                                                      -65.8980199769778, -60.4475595973147, -67.7543165093134, -60.4616894158005, 
                                                                      -67.9720957068844, -62.2184680275876, -66.2345680342004, -64.1523442367459, 
                                                                      -62.5435163857161, -65.9127866479611, -66.8874734854608, -62.0859917484373, 
                                                                      -66.8762861503705), Catcht = c(18L, 95L, 32L, 40L, 16L, 49L, 
                                                                                                     22L, 14L, 86L, 45L, 3L, 51L, 45L, 41L, 19L)), row.names = c(NA, 
                                                                                                                                                                 -15L), class = "data.frame")


# Start the analysis
an <- getData("GADM", country = "ATA", level = 0)
e <- extent(-70,-40,-68,-60)
rc <- crop(an, e)
proj4string(rc) <- CRS("+init=epsg:4326")

rc3 <- st_as_sf(rc) 

# Don't think you need st_shift_longitude, removed
catk <- st_as_sf(catk, coords = c("Longitude", "Latitude"),  crs = 4326) 


Grid <- rc3 %>%
  st_make_grid(cellsize = c(1,0.4)) %>% # para que quede cuadrada
  st_cast("MULTIPOLYGON") %>%
  st_sf() %>%
  mutate(cellid = row_number()) 




result <- Grid %>%
  st_join(catk) %>% 
  group_by(cellid) %>% 
  summarise(sum_cat = sum(Catcht))



ggplot() +
  geom_sf(data = Grid, color="#d9d9d9", fill=NA) +
  # Add here results with my mock data by grid
  geom_sf(data = result %>% filter(!is.na(sum_cat)), aes(fill=sum_cat)) +
  geom_sf(data = rc3) +
  theme_bw() + 
  coord_sf() +
  scale_alpha(guide="none")+
  xlab(expression(paste(Longitude^o,~'O'))) +
  ylab(expression(paste(Latitude^o,~'S')))+
  guides( colour = guide_legend()) +
  theme(panel.background = element_rect(fill = "#f7fbff"),
        panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank())+
  theme(legend.position = "right")+
  xlim(-69,-45)

^{Created on 2022-03-23 by the reprex package (v2.0.1)}

Yeah!!! you solve this!!!. And yes, my trouble is with `st_shift_longitude()`. I was follow another example and keep this function. Thanks @dieghernan !!! — MauroMardones, Mar 24 '22 at 14:42

How can plot my own data in a grid in a map sf but return vacum

1 Answers1