BookmarkSubscribeRSS Feed
🔒 This topic is solved and locked. Need further help from the community? Please sign in and ask a new question.
jonathanch
Calcite | Level 5

Good afternoon,

I'd like to know how could we use proc GMAP after proc GAM ?

I've looked this ressource:

http://www.sas.com/offices/NA/canada/en/resources/asset/McGahan-Maps.pdf

I didn't manage to execute a corectly program...

Thanks a lot for your help.

1 ACCEPTED SOLUTION

Accepted Solutions
PGStats
Opal | Level 21

L'ensemble de vos données (contrats) peut servir à l'ajustement du modèle GAM. Vous spécifiez cet ensemble de données dans PROC GAM DATA=mesDonnees;

Vous assemblez aussi un autre ensemble de données décrivant des clients fictifs de 40 ans situés au centroïde de chaque département. Vous réclamez les prédictions de votre modèle pour ces clients fictifs dans SCORE DATA=mesClientsFictifs OUT=mesPredictions; Ce sont ces prédictions que vous illustrez sur votre carte.

En résumé,

le dataset mesDonnees contiendrait les variables X, Y, ageconducteur, Nombre

et le dataset mesClientsFictifs contiendrait les variables ID_GEOFLA, X(=X_CENTROID), Y(=Y_CENTROID), ageconducteur(=40)

La procédure créera le dataset mesPredictions contenant les variables ID_GEOFLA, X, Y, ageconducteur(=40), P_Nombre

où seules les variables ID_GEOFLA et P_Nombre sont requises par GMAP pour colorer la carte.

proc gmap map=departements data=mesPredictions;

id ID_GEOFLA;

choro P_Nombre;

run;

quit;


Bonne chance!


PG (pas pédagogue pour deux sous Smiley Happy)

PG

View solution in original post

21 REPLIES 21
PGStats
Opal | Level 21

Please post the code you tried with some test data. - PG

PG
jonathanch
Calcite | Level 5

Ok thanks for advice me.

I'm working on a insurance data. This proc gam is an example with partial variables.

I have:

-geospatial variables such as longitud,lattitud, postal code...

-comportemental variables

-frequency of claims

-cost of claims

-exposure

I downloaded and imported the file FRANCE.shp.

I don't know how can I integrate my geospatial variables in PROC GAM...

I don't manage to make the link between PROC GAM and PROC GMAP.

My code:

%let geospatial=longitud lattitud postal_code;

PROC GAM DATA = ins.test;

MODEL nombre = spline(ageconducteur) &geospatial / dist = Poisson;

OUTPUT OUT=gam p uclm lclm PREDICTED;

RUN; QUIT;

PROC SORT DATA = gam NODUPKEY; BY ageconducteur; RUN; QUIT;

proc print data=gam;run;

DATA gam;

SET gam;

pred_nombre_gam = exp(P_nombre);

KEEP ageconducteur pred_nombre_gam &geospatial;

RUN;

*Import shp map;

proc mapimport out=commune datafile="C:\Users\roger\Documents\PROJET MJG\carto\GEOFLA\COMMUNES\commune.shp";

run;

data commune;

set commune;

cv=strip(CODE_DEPT!!CODE_CANT);

arr=strip(CODE_DEPT!!CODE_ARR);

rename CODE_DEPT=ID INSEE_COM=CODGEO X_CENTROID=X Y_CENTROID=Y;

*annotate;

DATA departement_anno ;

SET dep;

RETAIN xsys ysys "2" size 4 color "white" text "" function='label' when "after";

BY ID SEGMENT ;

IF FIRST.SEGMENT then function = "MOVE" ;

ELSE function = "DRAW" ;

keep id xsys ysys size color function x y text when;

RUN ;

run;

*proc gmap;

proc gmap data=??? map=commune_epci  annotate=departement_anno GOUT=MAPCHART ;

id CODGEO;

choro ???/ WOUTLINE=1 levels=5;

run;

PGStats
Opal | Level 21

What do you want the choropleth to look like?

To display the predicted values from your model for each dept, say, you would need to create a scoring dataset containing all your model predictors, including the polygon centroids, for each dept and use the SCORE statement in GAM to get the predictions. The dept code could then be used to make the link between your predictions and the map (assuming dept is also present in the shapefile). The centroids would only be used by GAM, not in GMAP.

This is a bit sketchy, I know. I don't have much experience with SAS mapping.

BTW, I don't understand why you would use option NODUPKEY when sorting predictions. That will keep only one record from each age group.

PG

PG
jonathanch
Calcite | Level 5

Good afternoon,

It's not easy for me...I don't manage in step 2.

Step 3 and 4, correct me if i'm wrong, I red all your sketchy Smiley Happy.

A model for each department.

You can see the file .SHP and the table if you want.

https://wxs-telechargement.ign.fr/oikr5jryiph0iwhw36053ptm/telechargement/inspire/GEOFLA_THEME-DEPAR...

Or:

My CENTROID1 and CONTRATS (insurance contracts) output:

                                        X          Y               SEGMENT       CODE_CHF     CODE_DEPT     NAME_DEPT

1

919315

6541572

1

053

01

AIN

2

918646

6540637

1

053

01

AIN

3

919271

6538851

1

053

01

AIN

4

919073

6537700

1

053

01

AIN

5

918313

6534729

1

053

01

AIN

6

917787

6532870

1

053

01

AIN

7

917348

6531205

1

053

01

AIN

8

917003

6529415

1

053

01

AIN

9

916337

6528661

1

053

01

AIN

10

916257

6525769

1

053

01

AIN

          ncontrat    exposition zone ... ...      ageconducteur                                                                                   freq                                       CODE_DEPT

12170.74A5431643D2181RC0.001-0.30111137
27090.18B78221002E2601RC0.001-1.71480278
37140.48C90326112E41134BG687.821-0.73397394
48520.27F753910012E1102DO96.641-1.30933486
510830.51E40495012E31132DO70.881-0.67334579

freq=number of claims  (dependant variable)

Step 1: I import the file SHP:

proc mapimport out=departement datafile="C:\Users\roger\Documents\PROJET MJG\carto\GEOFLA\DEPARTEMENTS/departement.shp";

run;

Step 2: many to many MERGE by code_dep to get my  CONTRATS table with X and Y (longitud,lattitude).

There are eight coordonate (I suppose town hall,prefecture...) by district for example so I don't manage. Proc sql but how? I find always many missing values.

Step 3: PROC GAM thanks to your advice.

PROC GAM DATA = merge_table;

MODEL freq = spline(ageconducteur)   spline(Y) spline(x) / dist = Poisson;

SCORE data=centroid1 out=estimate ;

BY code_dep;

RUN; QUIT;

Step 4: proc gmap

I drop X and Y (centroids) then proc gmap:

proc gmap data=centroid1 map=departement  annotate=departement_anno GOUT=MAPCHART ;

id code_dept;

choro myestimators_GAM / WOUTLINE=1 levels=5;

run;

If I don't respect the clauses on this website, I'm sorry to the moderator.

Thanks again PG for your help.

PGStats
Opal | Level 21

Poursuivons en français, si vous le voulez bien.

L'énoncé BY code_dep; dans l'appel de proc GAM réclame l'ajustement d'un modèle distinct pour chaque département. Ce n'est sans doute pas ce que vous souhaitez.

Vos données de fréquence des réclamations (et d'âge moyen des conducteurs) ont-elles une résolution spatiale plus fine que le département? Si oui, leurs coordonnées sont-elles compatibles avec les X-Y du fichier SHP?

Pour prédire des fréquences il faudra spécifier tous les prédicteurs. Pour l'âge par exemple, vous avez plusieurs choix. Vous pourriez utiliser l'âge moyen des conducteurs ou l'âge moyen de vos clients ou produire une carte pour chaque tranche d'âge, etc. Vous pourriez vouloir une carte distincte pour les hommes et les femmes. Ces choix détermineront le contenu du jeu de données que vous mentionnerez dans la commande SCORE.

PG

PG
jonathanch
Calcite | Level 5

Ok je tiens compte des différents choix possibles. Pour l'instant je m'entraîne sur une table dont les predicteurs ont des spécificités communes. Comme vous le voyez, je n'est pas vraiment d'experience SAS/MAP.

Cependant oui je peux avoir une résolution spatiale plus fine grâce au fichier commune (commune.SHP)

Ok pour ne pas utiliser BY code_dept pour avoir un seul et même modèle.

Là où je suis bloqué c'est pour affecter X et Y à ma table principale CONTRATS. Les X et Y étant dans le fichier SHP.

Si vous regardez l'extrait des tables (ci-dessus), vous verrez qu'il y a de nombreuses coordonnées par département d'où mon problème.

Avez-vous une idée?

PGStats
Opal | Level 21

Le plus précis, à priori, serait de géocoder les adresses à partir du code postal. Certaines ressources sont disponibles à cette fin (voir Europe - Postcode Point Lists, Zipcodes and Geocoders with Type: Postcode Geocoders and Country: Fra...) ainsi que proc geocode. Mais cela est probablement trop précis pour vos besoins.

Autrement, vous pouvez lire la table de données (DBF) qui accompagne le fichier SHP. Par exemple

filename dbDept "&sasforum\datasets\DEPARTEMENT.dbf";

proc dbf db5=dbDept out=deptData; run;

proc print data=deptData(obs=10); id id_GEOFLA; run;

et associer les X_CENTROID et Y_CENTROID qu'elle contient pour chaque département à vos adresses. Le problème supplémentaire que cela pose par rapport aux codes postaux est l'orthographe des noms (de département, de commune, etc) qui peut, comme vous le savez, être très variable. 

PG

PG
jonathanch
Calcite | Level 5

Bonjour,

Voilà le dernier résultat que j'obtiens après un proc GAM:

PROC GAM DATA = affectation;

MODEL nombre = spline(ageconducteur) spline(x) spline(y) / dist = Poisson;

OUTPUT OUT=gam p uclm lclm PREDICTED;

score data=commune_ageco out=estimate;

id code_com;

RUN; QUIT;

proc gmap data=estimate map=commune   GOUT=MAPCHART ;

id code_com;

choro P_ageconducteur / levels=4;

run;

Extrait table ESTIMATE:

     code_com    X                   Y                     P_nombre      P_ageconducteur  P_x                P_y                           LINP_nombre

1100037547678883200.15024-0.047510.033764-0.015031-1.89552
2100178045682073800.14884-0.042570.016636-0.015217-1.90488

gmap26.png

Je n'ai pas une map entière mais uniquement des fragments. Une erreur de ma part ou bien la map de GEOFLA qui n'est pas complète.

Pourtant un professeur a reussi avec les map de GEOFLA: https://www.google.fr/url?sa=t&rct=j&q=&esrc=s&source=web&cd=2&cad=rja&uact=8&ved=0CDEQFjAB&url=http...

Beaucoup de possibilité avec proc gmap mais je pense ne pas tout comprendre, pourtant j'aurais essayé.

Acheter environ 600 euros une carte c'est pas possible Smiley Happy

PGStats
Opal | Level 21

Bonjour

C'est un bon début! L'allure de la carte donne l'impression qu'elle ne montre que les communes pour lesquelles vous avez des contrats. De quelle façon la table commune_ageco a-t-elle été créée?

J'arrive à produire une carte complète (des départements, illustrant des données aléatoires) à l'aide du programme suivant:

/* Folder &sasforum\datasets contains files DEPARTEMENT.DBF, DEPARTEMENT.SHP and

DEPARTEMENT.SHX */

proc mapimport out=depts datafile="&sasforum\datasets\DEPARTEMENT.SHP"; run;

/* Associate a predicted number (p_nombre) with each Departement */

data test;

set depts; by code_dept;

if first.code_dept;

p_freq = ranpoi(123456, 10);

run;

/* Draw a choropleth of the predicted numbers */

proc gmap map=depts data=test;

id code_dept;

choro p_freq;

run;

quit;

PG

PG
jonathanch
Calcite | Level 5

Encore plus intéressant SAS XML MAPPING...mais j'y connais rien!

jonathanch
Calcite | Level 5

Ah oui les coordonnées de la carte sont complètes...

La table commune_ageco contient les X/Y et la variable ageconducteur pour faire fonctionner PROC GAM. C'est cette table qui me pose des problème car j'ai des valeurs manquantes dans cette table..

jonathanch
Calcite | Level 5

Ok vous prenez la première observation à chaque département dans le fichier departement.shp. Je vais essayer et vous tiens au courant. Thanks!

PGStats
Opal | Level 21

Pour le SCORE de proc GAM vous devez produire un dataset de tous les départements (ou communes) avec les X, Y et ageconducteur pour lesquels vous souhaitez obtenir une valeur prédite. Les valeurs manquantes dans les prédicteurs entraîneront des valeurs manquantes dans les valeurs prédites. Je vous conseille d'utiliser le fichier DBF comme base pour constituer le dataset DATA= de l'énoncé SCORE car il vous donne la position des centroides, alors que la première observation d'un fichier SHP ne vous donne que la position d'un point sur la frontière.

PG
jonathanch
Calcite | Level 5

gmap.png

Bien mieux Smiley Happy

Par contre quand je n'utilise pas l'option STATISTIC=freq... j'obtiens "WARNING: Some observations were discarded when charting P_nombre.  Only first matching observation was used.". Est-ce très important? Smiley Happy

Je voudrais aussi pouvoir déterminer de façon moins arbitraire le nombre de "pattern" (bornes de couleur et leurs nombres). Je me suis dit pourquoi pas utiliser une classification automatique ou encore kmeans method pour extraire le nombre de clusters qui seront bornés..

Cela serait une bonne idée?

Voilà dernière petite question car je vous ai demandé déjà beaucoup. Bonne continuation à vous.

sas-innovate-2024.png

Join us for SAS Innovate April 16-19 at the Aria in Las Vegas. Bring the team and save big with our group pricing for a limited time only.

Pre-conference courses and tutorials are filling up fast and are always a sellout. Register today to reserve your seat.

 

Register now!

What is ANOVA?

ANOVA, or Analysis Of Variance, is used to compare the averages or means of two or more populations to better understand how they differ. Watch this tutorial for more.

Find more tutorials on the SAS Users YouTube channel.

Discussion stats
  • 21 replies
  • 2888 views
  • 6 likes
  • 2 in conversation