Seismic Data – fixed?

As described in my last post, I was seeing significant gaps in the seismic event data I was retrieving from the INGV service. So I re-read their docs. Silly me, I’d missed the option to include query arguments restricting the geo area of the events (I had code in a subsequent script doing this).

While tweaking the code to cover these parameters I also spotted a really clumsy mistake. I had a function doing more or less this –

for each event element in XML DOM:
        extract event data
        add event data to list
        return list

D’oh! Should have been –

for each event element in XML DOM:
        extract event data
        add event data to list
return list

I’ve also improved error handling considerably, discriminating between genuine HTTP errors and HTTP 204 No Content. Now I’ve narrowed the geo area and reduced the time window for each GET down to 1 hour, there are quite a lot of 204s.

I’m now running it over the time period around the l’Aquila quakes as a sanity check. Jeez, 20+ events in some hours, 10+ in most.

Assuming this works ok, I’ll run it over the whole 1997-2017 preiod, hopefully in ~12 hours time I’ll have some usable data.

PS. Looking good, for the 30 days following that of the l’Aquila big one, it produced:

in_zone_count = 8877
max_depth = 62800.0
max_magnitude = 6.1

 

 

 

Advertisements

Seismic Data Wrangling

Following my interim plan of proceeding software-only (until I’ve the funds to get back to playing with hardware), I’ve been looking at getting seismic event data from the INGV Web Service into a Keras/Tensorflow implementation of PredNet.

My code is on GitHub, and rather than linking to individual files which I may rename, I’ll put a README over there with pointers.

As a first step, I put together code to pull of the data and dump it down to simple CSV files. This appeared to be working. The demo implementation of PredNet takes HDF5 data from the KITTI  vision dataset (videos from a car on road around Karlsruhe), extracting it into numpy arrays, with the PredNet engine using Keras. To keep things simple I wanted to follow the same approach. I’m totally new to HDF5 so pinged Bill Lotter of the PredNet project for clues. He kindly gave me some helpful tips, and concurred with what I’d been thinking – keep the CSV data, process that into something PredNet can consume.

The data offered by the Web Service is good XML delivered over HTTP (props to INGV). But it does include a lot of material (provenance, measurement accuracy etc) that isn’t needed here. So my service-to-CSV code parses out just the relevant parts, producing a line for each event:

datetime, latitude, longitude, depth, magnitude

e.g.

2007-01-02T05:28:38.870000, 43.612, 12.493, 7700, 1.7

I couldn’t find the info anywhere, but it appears that the INGV service records go back at least to somewhere in the early 1990’s, so I chose 1997-01-01T00:00:00 as a convenient start datetime, giving me 20 years of events.

For this to be a similar shape to what PredNet expects, I will aggregate events within a particular time period (actually I think taking the most significant event in that period). I reckon 6 hour periods should be about right. This also seemed a reasonable window for calling the service (not). I’ll filter down the events to just those within the region of interest (northern Italy, see earlier post)  and then scale the latitude & longitude to an easy integer range (probably [128, 128]). For a first pass I’ll ignore the depth field.

As it happens, I’m well on the way to having implemented this. But along the way I did a few sanity checks, eg. checking for maximum event magnitude in the region of interest, (I got 4.1), and it turned out I was missing some rather significant data points.  Here’s one I checked for:

The 2009 L’Aquila earthquake occurred in the region of Abruzzo, in central Italy. The main shock occurred at 03:32 CEST (01:32 UTC) on 6 April 2009, and was rated 5.8 or 5.9 on the Richter magnitude scale and 6.3 on the moment magnitude scale; its epicentre was near L’Aquila, the capital of Abruzzo, which together with surrounding villages suffered most damage.

Nope, it wasn’t in the CSV, but the Web Service knows all about it:

http://webservices.ingv.it/fdsnws/event/1/query?eventId=1895389

Doing a service call for that whole day:

http://webservices.ingv.it/fdsnws/event/1/query?starttime=2009-04-06T00:00:00&endtime=2009-04-06T23:59:59

–  yields 877 events – nightmare day!

I’d set the timeout on the HTTP calls to 2 seconds, but there is so much data associated with each event that this was woefully inadequate. Since upped to 5 mins.

Manually checking calls, I was also sometimes getting a HTTP Status Code of 413 Request Entity Too Large. This puzzled me mightily – still does actually. It says request entity, not requested (or response) entity, but the way it’s behaving is that the response requested is too large. Either way I reckon the spec (latest is RFC7231) is a little open to misinterpretation here. (What the heck – I’ve mailed the IEFT HTTP list about it – heh, well well, I’ve co-chaired something with the chair…).

Anyhow, I’ve also tweaked the code to make calls over just 1 hour windows, hopefully it’ll now get the stuff it was missing.

Hmm…I’ve got it running now and it’s giving errors throughout the year 2000, which should be trouble-free. I think I’ll have to have it make several passes/retries to ensure I get the maximum data available.

Drat! It’s giving me Entity Too Large with just 1 hour windows, e.g.

http://webservices.ingv.it/fdsnws/event/1/query?starttime=2000-12-13T01:00:00&endtime=2000-12-13T02:00:00

I need to fix this…

 

 

 

 

Candidate Neural Network Architecture : PredNet

While I sketched out a provisional idea of how I reckoned the network could look, I’m doing what I can to avoid reinventing the wheel. As it happens there’s a Deep Learning problem with implemented solutions that I believe is close enough to the earthquake prediction problem to make a good starting point : predicting the next frame(s) in a video. You train the network on a load of sample video data, then at runtime give it a short sequence and let it figure out what happens next.

This may seem a bit random, but I think I have good justification. The kind of videos people have been working with are things like human movement or motion of a car. (Well, I’ve seen one notable, fun, exception : Adversarial Video Generation is applied to the activities of Mrs. Pac-Man). In other words, a projection of objects obeying what is essentially Newtonian physics. Presumably seismic events follow the same kind of model. As mention in my last post, I’m currently planning on using online data that places seismic events on a map – providing the following: event time, latitude, longitude, depth and magnitude. The video prediction nets generally operate over time on x, y with R, G, B for colour. Quite a similar shape of data.

So I had a little trawl of what was out there.  There are a surprisingly wide variety of strategies, but one in particular caught my eye : PredNet. This is described in the paper Deep Predictive Coding Networks for Video Prediction and Unsupervised Learning (William Lotter, Gabriel Kreiman & David Cox from Harvard) and has supporting code etc. on GitHub. Several things about it appealed to me. It’s quite an elegant conceptual structure, which translates in practice into a mix of convnets/RNNs, not too far from what I anticipated needing for this application. This (from the paper) might give you an idea:

prednet-block

Another plus from my point of view was that the demo code is written using Keras on Tensorflow, exactly what I was intending to use.

Yesterday I had a go at getting it running.  Right away I hit a snag: I’ve got this laptop set up for Tensorflow etc. on Python 3, but the code uses hickle.py, which uses Python 2. I didn’t want to risk messing up my current setup (took ages to get working) so had a go at setting up a Docker container – Tensorflow has an image. Day-long story short, something wasn’t quite right. I suspect the issues I had related to nvidia-docker, needed to run on GPUs.

Earlier today I decided to have a look at what would be needed to get the PredNet code Python3-friendly. Running kitti-train.py (Kitti is the demo data set) led straight to an error in hickle.py. Nothing to lose, had a look. “Hickle is a HDF5 based clone of Pickle, with a twist. Instead of serializing to a pickle file, Hickle dumps to a HDF5 file.“. There is a note saying there’s Python3 support in progress, but the cause of the error turned out to be –

if isinstance(f, file):

file isn’t a thing in Python3. But kitti-train.py was only passing a filename to this, via data-utils.py, so I just commented out the lines associated with the isinstance. (I guess I should fix it properly, feed back to Hickle’s developer.)

It worked! Well, at least for kitti-train.py. I’ve got it running in the background as I type. This laptop only has a very wimpy GPU (GeForce 920M) and it took a couple of tweaks to prevent near-immediate out of memory errors:

export TF_CUDNN_WORKSPACE_LIMIT_IN_MB=100

kitty_train.py, line 35
batch_size = 2 #was 4

It’s taken about an hour to get to epoch 2/150, but I did renice Python way down so I could get on with other things.

Seismic Data

I’ve also spent a couple of hours on the (seismic) data-collecting code. I’d foolishly started coding around this using Javascript/node, simply because it was the last language I’d done anything similar with. I’ve got very close to having it gather & filter blocks of from the INGV service and dump to (csv) file. But I reckon I’ll just ditch that and recode it in Python, so I can dump to HDF5 directly – it does seem a popular format around the Deep Learning community.

Radio Data

Yes, that to think about too.

My gut feeling is that applying Deep Learning to the seismic data alone is likely to be somewhat useful for predictions. From what I’ve read, the current approaches being taken (in Italy at least) are effectively along these lines, leaning towards traditional statistical techniques. No doubt some folks are applying Deep Learning to the problem. But I’m hoping that bringing in radio precursors will make a major difference in prediction accuracy.

So far I have in mind generating spectrograms from the VLF/ELF signals. Which gives a series of images…sound familiar? However, I suspect that there won’t be quantitatively all that much information coming from this source (though qualitatively, I’m assuming vital).  As a provisional plan I’m thinking of pushing it through a few convnet/pooling layers to get the dimensionality way down, then adding that data as another input to the  PredNet.

Epoch 3/150 – woo-hoo!

PS.

It was taking way too long for my patience, so I changed the parameters a bit more:

nb_epoch = 50 # was 150
batch_size = 2 # was 4
samples_per_epoch = 250 # was 500
N_seq_val = 100 # number of sequences to use for validation

It took ~20 hours to train. For kitti_evaluate.py, it has produced some results, but also exited with an error code. Am a bit too tired to look into it now, but am very pleased to get a bunch of these:

 

 

 

New Strategy for Seismic Data

I’m a massive procrastinator and not a very quick thinker. There is something positive about all this. I’ve barely looked at code in the past few weeks, but have been thinking about it, and may well have saved some time. (I usually get there in the end…).

The neural net setup I had in mind was based on the assumption that I’d have my own, local, data sources (sensors). But hardware is still on hold until I find some funds. So I’ve been re-evaluating how best to use existing online data.

Now I am pleased with the idea of taking the VLF radio data as spectrograms and treating them (conceptually) as images, so I can exploit existing Deep Learning setups. If I’m not getting the seismic data as time series from a local sensor(s) but from INGV, I can use the same trick. They have a nice straightforward Web Service offering Open Data as QuakeML (XML) over HTTP from URL queries. They also render it like this:

seismo-map

So I’m thinking of taking the magnitude & depth data from the web service and placing it in a grid (say 256×256) derived from the geo coordinates. To handle the time aspect, for each cell in the grid I reckon picking the max magnitude event over each 6 (?) hour window. And then (conceptually) treating this as a pixel map. I need to read up a little more, but this looks again like something that might well yield to convnets, constructed very like the radio data input.

PS.

I made a little start on coding this up. First thing was to decide what area I wanted to monitor. Key considerations were: it must include here (self-preservation!); it must cover a fair distance around the VLF data monitoring station who’s data I’m going to use (Renato Romero’s, here); it should include the main seismically active regions likely to impact those two places.

You can get an idea of the active regions from the map above, here’s another one showing estimated risk:

hotspots

In one of the papers I’ve read that the radio precursors are only really significant for a 100km or so (to be confirmed), so I’ve chosen the area between the 4 outer markers here:

area

This corresponds to latitude 40N-47N, longitude 7E-15E.

The marker middle-left is the VLF station, the one in the middle is where I live.

It looked like the kind of thing where I could easily get my lats & longs mixed up so I coded up the little map using Google Maps, was very easy, source – rendered it on JSFiddle.

Next I need to get the code together to make the services requests, filter & aggregate the seismic event data is some convenient for ready for network training.

Morse Code Practice Key

Recently I’ve been trying to fill in some of the massive holes in my knowledge about radio. For this reason I’d quite like to have a go at the radio amateur license exams. No idea how to go about it though, the geographic complications. I suppose I might stand a chance with the Italian version, as long as I could take a dictionary with me. Anyone know anything about this?

Anyhow, even though it’s rather anachronistic (and no longer a requirement for the exams), the radio amateur sites all have some mention of Morse Code. I’ve always wanted to learn it, I guess from watching too many spy films. There are loads of things I should be doing, but yesterday’s procrastination was making this gadget. Good fun.

Hardware Delusions

The key front end sensors I wish to build for this project are an ELF/VLF radio receiver and a seismometer. The frequency ranges of interest here are < 20kHz,  in other words, in the audio range (and probably extending a little lower).

As it happens, in a past life I studied audio frequency electronics, transducers, signals and systems, DSP and so on formally for 3 years, have written articles on (nonlinear) analog electronics for a magazine, and probably more significantly have been an electronic music hobbyist for around 40 years. In short, I consider myself something of an expert in the field. The word for this is hubris.

I started planning the sensors like a bull at a gate. On the seismic side, I hadn’t really thought things through very well. On the radio side – I’d only really skimmed Radio Nature, and my knowledge of radio reception is minimal. Since then, the flaws in my ideas have poured out.

Seismic Errors

I’ve got a design for seismic signal sensors roughed out. While a magnet & coil is a more traditional way of detecting audio frequency deflections, I thought it would be neater somehow to use semiconductor Hall Effect devices. A standard design for a proximity detector is one of these components (which are housed much like transistors) backed by a magnet. When a something like a piece of iron passes by, the magnetic flux varies and hence the output of the device (linear output devices are available).

So for my seismometer, the moving part will be a steel ball bearing on a spring, hanging in a jar of oil (for damping). There will be 3 sensors located in the x, y & z directions (N-S, E-W, up-down) relative to this.

One potential complication with this setup had occurred to me. For a (relatively) linear response, the ball bearing would have to move in line with the face of the sensor. Obviously, in practice, most of the time the movement will be off-axis. However, my thinking went, there should still be enough information coming from all 3 sensors in combination to potential determine the deflection of the ball bearing. The data produced by these sensors will ultimately go into a neural network system, and they’re good at figuring out peculiar relationships.

But I’d missed another potential source of problems, it only came to me a couple of days ago. There is likely to be significant, pretty complex, interaction between the ball bearing and all 3 magnets. Whether or not this additional complication will be such that that the directional seismic information is totally obfuscated remains to be seen. I plan to experiment, maybe I’ll need 3 independent sensors…

Loopy

A little learning is a dang’rous thing. The danger I faced was wasting time & money in building a VLF receiver that simply couldn’t work.

I’d only skimmed the material, but something about the use of a coil as a receiver appealed to me. But the designs I’d seen were all pretty big, say around 1m in diameter. Hmm, thinks I, why not shrink that down to say 30cm and just boost the gain of the receiver circuit. It was only after I’d wound such a coil and picked up nothing but hum & noise that I got around to reading a little more.

It turns out there are two related issues involved: the way a small (relative to wavelength) loop antenna works isn’t exactly intuitive, and also its output is very low. It’s frequency-dependent, but the level of the desired signal is at a similar order of magnitude as the thermal noise generated by the loop, less than that of many op amps. The good Signore Romero, author of Radio Nature, has a practical discussion of this in his description of A Minimal ELF Loop Receiver. (Being at the low end of the frequency range of interest make this rather a worst-case scenario, but the points still apply). Basically there’s a good reason for having a big coil.

Another possible design flaw coming from my lack of learning is that I initially thought it would make sense to have coils in the x, y & z dimensions. As it turns out, because VLF signals are propagated as ground waves (between the surface of the planet and the ionosphere), pretty much all a coil in the horizontal plane will pick up is local noise such as mains hum. But I’m not yet discarding the inclusion of such a loop. Given the kind of neural net processing I have in mind, a signal that is comprised of little more than local noise may well be useful (in effect subtract this from the other signals).

But even having said all this, a loop antenna may still be of no use for me here – Noise Annoys. Renato has an image that nicely sums up the potential problem:

renato-noise

Right now I don’t have the funds to build a loop antenna of any description (big coils use a lot of wire!) but as and when I can, I’ll probably be looking at something along the lines of Renato et al’s IdealLoop (the image above comes from that page).

I do have the components to put together some kind of little portable whip antenna (electric field) receiver, I think I’ll have a look at that next, particularly to try and get an idea of how the noise levels vary in this locale.

I’ve also got one linear Hall effect sensor, so I can have a play around with that to try and get some idea of my seismometer design’s viability.

 

Summary for my Mother

She doesn’t have a science/technology/coding background, but asked me. But she does understand purple prose, mailed her this:

This is my view : it has been observed that many earthquakes have been preceded by unusual electrical activity in the atmosphere, and this has been made evident in the Extremely Low Frequency/Very LF radio bands.

These happen to be the frequency bands that audio works at, so a lot of existing technologies can be applied to its analysis (although the signals are radio rather than acoustic). There are plausible scientific explanations for these phenomena – if virtually any kind of rock is heavily compressed or stretched it generates electricity, at least some of which is expressed as radio waves. Rocks get squeezed and pulled before a split. Research to date has generally been inconclusive, even negative on this having potential for predicting earthquakes.

The biggest problem in my opinion is that there is a massive amount of unrelated noise in the system – things like lightning strikes are very loud in the wavebands in question, let alone solar flares. But recent developments in machine learning, aka Deep Learning, have been very good at pulling out features of interest from a very confused situation.

My hypothesis is that modern neural nets will be able to detect patterns well enough (beyond what traditional techniques can find) to provide an early warning system. Don’t get me wrong, I’m only half convinced myself, but if there’s a 10% of the idea making sense, 10% chance of it being codable (that bit’s 100%, but hypothetically) and a 10% chance of being able to exploit it, I’ll be out in the garden in a tent and will live to talk bollocks again. Or come back in shame-faced – but still alive.