/* SAS program to accompany "The path of zip codes" by Rick Wicklin published 18SEP2017 on The DO Lopp blog http://blogs.sas.com/content/iml/2017/09/18/path-of-zip-codes.html The program requires SASHELP.ZIP data set, which is distributed with SAS MAPSGFK.US_STATES data set, which I think is only distributed if you have a SAS/GRAPH licence */ /* Sort zip codes and exclude territories and special zip codes */ proc sort data=sashelp.zipcode(where=(StateCode NOT IN ("PR", "FM", "GU", "MH", "MP", "PW", "VI") /* exclude territories */ AND ZIP_Class = " ")) /* exclude "special" zip codes */ out=zipcode(keep=ZIP x y City StateCode); by zip; run; ods graphics / width=800px height=450px antialiasmax=50000 labelmax=10000; title "Path of US Zip Codes in Numerical Order"; proc sgplot data=zipcode noborder; where StateCode NOT IN ("AK", "HI"); /* contiguous US */ series x=X y=Y / group=StateCode; xaxis display=none; yaxis display=none; run; /* To see only a subset, use a WHERE clause, such as where StateCode IN ('UT', 'CO', 'AZ', 'NM'); */ /********************/ /* create a map of a subset of states (zoom in). The following code not only extracts the data, but uses some SAS/IML code to set the width/length of the resulting graph so that short wide states (TN, NC) and long tall states (CA, FL) are visualized without being overly stretched. */ data US; set mapsgfk.US_States(keep=segment state statecode long lat); where StateCode NOT IN ("AK", "HI", "PR"); by State Segment; if first.Segment then PolyID+1; BorderStateCode = StateCode; drop StateCode; run; %let myStates = 'IL', 'IA', 'MO'; %let myStates = 'UT', 'CO', 'AZ', 'NM'; data TMP; set zipcode(where=(StateCode in (&myStates))) US(where=(BorderStateCode in (&myStates))); run; /* compute YPIXELS for a 640 x YPIXEL image */ proc iml; use TMP; read all var {"Lat" "Long"}; close; XRange = range(Long); YRange = range(Lat); XPixels = 640; YPixels = round( XPixels * YRange / XRange ); print XPixels YPixels; call symputx("YPixels", YPixels); QUIT; ods graphics / width=640px height=&YPixels.px; title2 "StateCode in (&myStates)"; proc sgplot data=TMP noborder noautolegend; series x=X y=Y / group=StateCode; polygon x=Long y=Lat ID=PolyID / fill fillattrs=(transparency=0.75) outline label=BorderStateCode labelattrs=(Size=14 weight=bold); xaxis display=none; yaxis display=none; run; /********************/ /* create the Sectional Center Facility (SCF) value as '123xx' */ data Zip1; length SCF $5.; set zipcode; FloorZip = floor(zip/100); /* round down to nearest 100 */ SCF = putn(FloorZip,"Z3.")||"xx"; /* Sectional Center Facility, eg, 841xx */ keep x y zip StateCode City SCF; run; proc freq data=Zip1; where StateCode='UT'; tables SCF / nocum; run; /********************/ /* convenient macro to plot ONE state at a time and connect the SCFs (in order) within that state */ %macro ZipPath(StateCode, label=SCF); data State; set Zip1(where=(StateCode="&StateCode")); by SCF; if first.SCF then do; xx=x; yy=y; zipcount+1; path=zipcount; end; else do; xx=.; yy=.; path=.; end; run; data TMP; set State US(where=(BorderStateCode="&StateCode")); run; proc iml; use TMP; read all var {"Lat" "Long" "SCF"}; close; XRange = range(Long); YRange = range(Lat); XPixels = 640; YPixels = round( XPixels * YRange / XRange ); *print XPixels YPixels; call symputx("YPixels", YPixels); u = unique(SCF); Z = u[2] + " - " + u[ncol(u)]; call symputx("ZipRange", Z); QUIT; ods graphics / width=640px height=&YPixels.px; title "Path of Zip Codes for &StateCode"; title2 "&ZipRange"; proc sgplot data=TMP noborder noautolegend; format path 2.; scatter x=X y=Y / markerattrs=(symbol=circlefilled) group=SCF transparency=0.4; series x=xx y=yy / lineattrs=(color=black thickness=2 pattern=solid); scatter x=xx y=yy / jitter datalabel=&label datalabelattrs=(weight=bold size=12); polygon x=Long y=Lat ID=PolyID / fill fillattrs=(transparency=0.8) outline label=BorderStateCode labelattrs=(Size=14 weight=bold); xaxis display=none; yaxis display=none; run; %mend; %ZipPath(UT); /* by default, label by SCF codes */ %ZipPath(NM); %ZipPath(SC, label=path); /* or label by 1,2,3,.... */ %ZipPath(FL, label=path); /*****************************/ /* graph of US that connects SCFs in order. Get rid of labels and make other simplifications to make the image cleaner. */ data State; set Zip1; where StateCode NOT IN ("AK", "HI"); /* contiguous US */ by SCF; if first.SCF then do; xx=x; yy=y; zipcount+1; path=zipcount; end; else do; xx=.; yy=.; path=.; end; run; data TMP; set State US; run; ods graphics / width=800px height=450px antialiasmax=50000 labelmax=10000; title "Path of Sectional Center Facilities within Zip Codes"; proc sgplot data=TMP noborder noautolegend; polygon x=Long y=Lat ID=PolyID / outline; series x=xx y=yy; xaxis display=none; yaxis display=none; run;