意匠の時空間的分析

基礎な手法を用いて日本における登録意匠の時空間的分析などを探検する目的で、このブログによって、手法などについてシェアしたり、意見交換したりします。

CSRの帰無仮設を棄却する検定手順(1)背景と分析地区(データクリーニング後)

  1. 背景

    点過程データを分析する際、その最初のステップとして、CSR(Complete Spatial Randomness)という帰無仮設を置き、最近隣距離法若しくは法格法を用いて、その仮設を棄却する確率を算出する。*1*2ここでは、Rのspatstatというパッケージ*3を用いて、佐賀県と福岡県の意匠登録者の空間点過程データをBivand, Pebesma, Gómez-Rubio (2008)が説明するDiggle(2003) の例に従って、分析する。

     

  2. 分析地区

    spパッケージのSpatialPointsオブジェクトに対して、spatstatがpppオブジェクトを利用するため、spatstatの機能を利用するようにspSagaFukuokaの地区を抽出した後、pppオブジェクトに変換しなければならない。

      • 佐賀県と福岡県の全地区(pppSagaFukuoka)

        
        > ## Bivand, Pebesma, Gómez-Rubio (2008)pp156-158
        > library(maptools)
        > summary(spSagaFukuoka)
        Object of class SpatialPoints
        Coordinates:
                min       max
        x 129.81595 131.11875
        y  33.00453  33.95653
        Is projected: FALSE 
        proj4string : [+proj=longlat +ellps=WGS84]
        Number of points: 1255
        > ## convert to scaled rectangle
        > spSagaFukuokaScaled <- elide(spSagaFukuoka, scale=TRUE, unitsq = FALSE)
        > summary(spSagaFukuokaScaled)
        Object of class SpatialPoints
        Coordinates:
          min       max
        x   0 1.0000000
        y   0 0.7307327
        Is projected: NA 
        proj4string : [NA]
        Number of points: 1255
        > pppSagaFukuoka <- as( spSagaFukuokaScaled, "ppp")
        > summary(pppSagaFukuoka)
        Planar point pattern:  1255 points
        Average intensity 1255 points per square unit
        
        Coordinates are given to 8 decimal places
        
        Window: rectangle = [0, 1] x [0, 1] units
        Window area = 1 square unit
        

        f:id:ttcrossroads09:20181127173045p:plain

      • 福岡市―筑紫野市(pppFukuokaChikushino)

        
        > ## done in earlier session
        > llCRS <- CRS("+proj=longlat +ellps=WGS84")
        > ## coordinates
        > Sr1 = Polygon(cbind(
        + c(130.191354,130.626400,130.626400,130.191354, 130.191354),
        + c(33.718652,33.718652,33.423697,33.423697,33.718652)
        + ))
        > Srs1 = Polygons(list(Sr1), "bbFukuokaChikushino")
        > spbbFukuokaChikushino = SpatialPolygons(list(Srs1), 1:1, proj4string = llCRS)
        > summary(spbbFukuokaChikushino)
        Object of class SpatialPolygons
        Coordinates:
               min       max
        x 130.1914 130.62640
        y  33.4237  33.71865
        Is projected: FALSE 
        proj4string : [+proj=longlat +ellps=WGS84]
        > spFukuokaChikushino <- spSagaFukuoka[spbbFukuokaChikushino,]
        > summary(spFukuokaChikushino)
        Object of class SpatialPoints
        Coordinates:
                min       max
        x 130.19690 130.58931
        y  33.42709  33.71665
        Is projected: FALSE 
        proj4string : [+proj=longlat +ellps=WGS84]
        Number of points: 601
        >
        

        f:id:ttcrossroads09:20181127174343p:plain

        
        > ## convert to scaled rectangle
        > plot(spFukuokaChikushino, pch=1)
        > spFukuokaChikushinoScaled <- elide(spFukuokaChikushino, scale=TRUE, unitsq = FALSE)
        > summary(spFukuokaChikushinoScaled)
        Object of class SpatialPoints
        Coordinates:
          min       max
        x   0 1.0000000
        y   0 0.7379017
        Is projected: NA 
        proj4string : [NA]
        Number of points: 601
        > pppFukuokaChikushino <- as(spFukuokaChikushinoScaled, "ppp")
        > summary( pppFukuokaChikushino)
        Planar point pattern:  601 points
        Average intensity 814.4716 points per square unit
        
        Coordinates are given to 8 decimal places
        
        Window: rectangle = [0, 1] x [0, 0.7379017] units
        Window area = 0.737902 square units
        > png(file="pppFukuokaChikushino.png", width=600, height=480)
        > plot( pppFukuokaChikushino, pch=1)
        > dev.off()
           
        

        f:id:ttcrossroads09:20181127174411p:plain

    • 福岡市の中央区博多区(pppChuoHakata)

      
      > ########## ChuoHakata 
      > Sr2 = Polygon(cbind(
      + c(130.354992,130.486107,130.486107,130.354992, 130.354992),
      + c(33.613987,33.613987,33.535556,33.535556,33.613987)
      + ))
      > Srs2 = Polygons(list(Sr2), "bbChuoHakata")
      > spbbChuoHakata = SpatialPolygons(list(Srs2), 1:1, proj4string = llCRS)
      > summary(spbbChuoHakata)
      Object of class SpatialPolygons
      Coordinates:
              min       max
      x 130.35499 130.48611
      y  33.53556  33.61399
      Is projected: FALSE 
      proj4string : [+proj=longlat +ellps=WGS84]
      > spChuoHakata <- spSagaFukuoka[spbbChuoHakata,]
      > summary(spChuoHakata)
      Object of class SpatialPoints
      Coordinates:
              min       max
      x 130.35518 130.48578
      y  33.53611  33.61334
      Is projected: FALSE 
      proj4string : [+proj=longlat +ellps=WGS84]
      Number of points: 348
      >
      

      f:id:ttcrossroads09:20181127174524p:plain

      
      > ## convert to scaled rectangle
      > spChuoHakataScaled <- elide(spChuoHakata, scale=TRUE, unitsq = FALSE)
      > summary(spChuoHakataScaled)
      Object of class SpatialPoints
      Coordinates:
        min       max
      x   0 1.0000000
      y   0 0.5913476
      Is projected: NA 
      proj4string : [NA]
      Number of points: 348
      > pppChuoHakata <- as(spChuoHakataScaled, "ppp")
      > summary(pppChuoHakata)
      Planar point pattern:  348 points
      Average intensity 588.4863 points per square unit
      
      Coordinates are given to 8 decimal places
      
      Window: rectangle = [0, 1] x [0, 0.5913476] units
      Window area = 0.591348 square units
      

      f:id:ttcrossroads09:20181127174558p:plain

    • 春日市(pppKasuga)

      
      > ## Earlier Saga-Fukuoka, Fukuoka-Chikushino Cities, Chuo-Hakata Wards 
      > ## subsets show high degrees of aggregation. But this is to be expected 
      > ## from terrain and population. This subset is from an area of Kasuga City  
      > ## Chikushino cities that has had a lot of planned development with 
      > ## that are newer developments with large areas of single family dwellings.
      > ## 33.556279, 130.481984
      > ## 33.516162, 130.409705
      > latNE = 33.556279
      > longNE = 130.481984
      > latSW = 33.516162
      > longSW = 130.409705
      > Sr1 = Polygon(cbind(
      + c(longNE ,longNE ,longSW ,longSW , longNE ),
      + c(latNE,latSW ,latSW ,latNE,latNE)
      + ))
      > Srs1 = Polygons(list(Sr1), "bbKasuga")
      > spbbKasuga = SpatialPolygons(list(Srs1), 1:1, proj4string = llCRS)
      > summary(spbbKasuga)
      Object of class SpatialPolygons
      Coordinates:
              min       max
      x 130.40971 130.48198
      y  33.51616  33.55628
      Is projected: FALSE 
      proj4string : [+proj=longlat +ellps=WGS84]
      > spKasuga <- spSagaFukuoka[spbbKasuga,]
      > summary(spKasuga)
      Object of class SpatialPoints
      Coordinates:
              min      max
      x 130.40993 130.4798
      y  33.51632  33.5561
      Is projected: FALSE 
      proj4string : [+proj=longlat +ellps=WGS84]
      Number of points: 50
      >
      

      f:id:ttcrossroads09:20181127174722p:plain

      
      > ## convert to scaled rectangle
      > spKasugaScaled <- elide(spKasuga, scale=TRUE, unitsq = FALSE)
      > summary(spKasugaScaled)
      Object of class SpatialPoints
      Coordinates:
        min       max
      x   0 1.0000000
      y   0 0.5692616
      Is projected: NA 
      proj4string : [NA]
      Number of points: 50
      > 
      > pppKasuga <- as(spKasugaScaled, "ppp")
      > summary(pppKasuga)
      Planar point pattern:  50 points
      Average intensity 87.83308 points per square unit
      
      Coordinates are given to 8 decimal places
      
      Window: rectangle = [0, 1] x [0, 0.5692616] units
      Window area = 0.569262 square units  
      

      f:id:ttcrossroads09:20181127174801p:plain