[R-sig-Geo] Shapefile and Basemap
sownal chand
@own@|ch@nd|m@ @end|ng |rom gm@||@com
Sun May 15 18:54:42 CEST 2022
Hi Micha
Thanks a lot for your assistance.
Have a blessed day
Sownalc
On Mon, May 16, 2022, 03:15 Micha Silver <tsvibar using gmail.com> wrote:
>
> On 5/15/22 03:52, sownal chand wrote:
> > Hi Micha,
> >
> > Thanks for your email and really appreciate your assistance.
> >
> > I see that there are alot of old versions of information on the net
> > with spatial analysis using R. I would like to request if you have any
> > particular site or reading materials for spatial analysis and
> > modelling in regards to R programming would really help.
> >
>
> One of the best online sources (in my opinion) is the book by Robin
> Lovelace "Geocomputation with R":
>
> https://geocompr.robinlovelace.net/
>
>
> Here's a good explanation of all the changes in the R-spatial packages,
> and why those changes were necessary:
>
> https://r-spatial.org/r/2022/04/12/evolution.html
>
>
> If you're following tutorials that depend on any of the packages: `sp`,
> `raster`, `rgdal`, then you should be aware that these will soon be
> deprecated and your procedure will not work in a year or so.
>
> (Authors of those tutorials will have to revise them...)
>
>
> Good luck
>
>
>
> > Once again, thanking you in advance
> > Kind regards
> > Sownalc
> >
> > On Sun, May 15, 2022, 03:50 Micha Silver <tsvibar using gmail.com> wrote:
> >
> > Hello
> >
> >
> > On 5/13/22 12:47, sownal chand wrote:
> > > Hello sir/madam,
> > >
> > > I am working with shape file of my country and the issue I am
> > facing
> > > is the shapefile is scattered while plotting it using basemap. I am
> > > using my sample point data which is attached to this email. I hope
> > > that some expert in this area would help in correcting the codes
> > below
> > > to show the shapefile in one location ( its a pacific centered map)
> >
> >
> > When you raised this question a few weeks ago, it was suggested to
> > avoid
> > the `sp` package with its SPDF data type, and instead focus on the
> > newer
> > `sf` package.
> >
> > There is also a replacement for the raster::getData function in
> > the new
> > `geodata` package.
> >
> >
> > Here is a much simpler version of what (I think) you are trying to
> > achieve:
> >
> >
> > # Load only three libraries
> >
> > library(sf)
> > library(tmap)
> > library(geodata)
> >
> >
> > # Read your list of data (You should remove the summary line in
> > advance...)
> >
> > dt <- read.csv("DataR.csv")
> > dt <- dt[complete.cases(dt),]
> > dt_sf <- st_as_sf(dt, coords=c("long", "lat"),
> > crs="EPSG:4326")
> > str(dt_sf)
> >
> > # Get Fiji boundary from geodata package
> > fiji <- gadm(country="FIJI", level=2, path=tempdir())
> > # Convert to sf object for tmap plotting
> > fiji <- st_as_sf(fiji)
> >
> >
> > # Visualize with tmap
> >
> > tmap_mode("view")
> > tm_basemap("OpenStreetMap.Mapnik") +
> > tm_shape(fiji) +
> > tm_borders(col="brown", lwd=2) +
> > tm_shape(dt_sf) +
> >
> > # size of symbols by yearly data. You can choose any
> > year, of
> > course
> >
> > tm_symbols(col="blue", size="Year2", scale=0.1)
> >
> >
> >
> > If you have a specific problem with a certain shapefile, you'll
> > have to
> > supply it to the list in order to get further help.
> >
> >
> > HTH
> >
> >
> > >
> > >
> >
> **************************************************************************************************
> > >
> > >
> > > library(sp)
> > > library(raster)
> > > library(rgdal)
> > > library(leaflet)
> > >
> > > #read.csv
> > > read.csv ("C://Users/Documents/data.csv") -> data.df
> > > head(data.df)
> > >
> > > hist(data.df$Year, breaks=20)
> > >
> > > #remove NA valuues in the spatial Data Frame
> > > data.df <- na.omit(data.df)
> > > View(data.df)
> > >
> > > plot(data.df$long, data.df$lat,
> > > ylab = "Latitude", xlab="Longitude") #boring!
> > >
> > > # Use the cex function to plot circle size as a function of a
> > variable
> > > plot(data.df$long, data.df$lat,
> > > cex = data.df$Year.7 * 0.045,
> > > ylab = "Latitude", xlab="Longitude")
> > >
> > > data.df_SPDF <- SpatialPointsDataFrame(coords = data.df[,c("long",
> > > "lat")],
> > > data =
> > > data.df[,c("Year", "Year.1", "Year.2","Year.3","Year.4")],
> > > proj4string =
> > > CRS("+init=epsg:4326")) # sets the projection to WGS 1984 using
> > > lat/long. Optional but good to specify
> > >
> > > # Summary of object
> > > data.df_SPDF
> > >
> > > # SPDFs partition data elements, e.g. the coordinates are stored
> > > separately from the data
> > > head(data.df_SPDF using coords)
> > >
> > > head(data.df_SPDF using data)
> > >
> > >
> > > # You can quickly access the data frame as per a standard data
> > frame, e.g.
> > > head(data.df_SPDF$Year)
> > >
> > >
> > > # You can use the plot or spplot function to get quick plots
> > > plot(data.df_SPDF)
> > >
> > > spplot(data.df_SPDF, zcol = "Year")
> > >
> > > FIJ_Adm_2 <- readOGR("FIJ_Adm_2_shapefile", "FIJ_Adm_2")
> > >
> > > # You first need the ISO3 codes for the country of interest.
> > >
> > > # The getData function then allows you to retrieve the relevant
> > admin
> > > level boundaries from GADM.
> > > FJI_Adm_2 <- raster::getData("GADM", country="FJI", level=2)
> > >
> > > # Plot both country and data points
> > > plot(FJI_Adm_2)
> > > points(data.df$long, data.df$lat,
> > > cex = data.df$Year * 0.045,
> > > ylab = "Latitude", xlab="Longitude",
> > > col="red")
> > >
> > > basemap <- leaflet() %>% addProviderTiles("CartoDB.Positron")
> > >
> > > basemap %>% addPolygons(data=FJI_Adm_2)
> > >
> > > # to change the colors/line weight
> > > basemap %>% addPolygons(data=FJI_Adm_2, color = "red",
> > > weight = 1, fillOpacity = 0.2)
> > >
> > >
> > > # If you want to add points as well
> > > basemap %>% addPolygons(data=FJI_Adm_2, weight = 2,
> > > popup = FJI_Adm_2$NAME_2) %>%
> > >
> > > addCircleMarkers(data=data.df_SPDF,
> > > color="red", radius = 2)
> > >
> > > library(wesanderson) # for a nice color palette
> > > colorPal <- colorNumeric(wes_palette("Zissou1")[1:5],
> > > data.df_SPDF$Year, n = 5)
> > >
> > >
> > > # colorPal is now a function you can apply to get the corresponding
> > > # color for a value
> > > colorPal(0.6)
> > >
> > > basemap %>% addPolygons(data=FJI_Adm_1, weight = 2, fillOpacity=0,
> > > popup = FJI_Adm_1$NAME_1) %>%
> > >
> > > addCircleMarkers(data=data.df_SPDF,
> > > color = colorPal(data.df$Year),
> > > radius = 2,
> > > popup = as.character(data.df$Year))%>%
> > >
> > > addLegend(pal = colorPal,
> > > title = "Average Temp for Year",
> > > values = data.df_SPDF$Year)
> > >
> > >
> > >
> > >
> >
> *************************************************************************************************
> > >
> > > Thanking you in advance
> > > sownalc
> >
> > --
> > Micha Silver
> > Ben Gurion Univ.
> > Sde Boker, Remote Sensing Lab
> > cell: +972-523-665918
> > https://orcid.org/0000-0002-1128-1325
> >
> --
> Micha Silver
> Ben Gurion Univ.
> Sde Boker, Remote Sensing Lab
> cell: +972-523-665918
> https://orcid.org/0000-0002-1128-1325
>
>
[[alternative HTML version deleted]]
More information about the R-sig-Geo
mailing list