Tutorial lecture on machine-learning methods for astronomers (classification, regression, feature engineering) at Astro Hack Week 2016, hosted at UC Berkeley and recorded by BIDS.
Exact day within Aug 29-Sep 2 approximate.
“Start thinking about machine learning as just another tool in your toolkit for doing inference… There's a time and a place for applying different tools.” – Joshua Bloom
“There's been this almost Cambrian explosion of different algorithms and approaches… it can often be very noisy and very difficult to figure out what's hype and what's actually very useful.” – Joshua Bloom
“You have to look at your data before you start throwing it into these frameworks. Just because these frameworks are nice and easy to use doesn't mean that you're allowed to get away with not actually looking at the data.” – Joshua Bloom
“If you get a perfect classifier and you see a confusion matrix that has no power off the diagonal and all the power in the diagonal, you've overfit. I guarantee it.” – Joshua Bloom
“There are lots of ways to hammer your thumb with your new tool. And you may think you're building this amazing house, but in the end, you've built complete crap… Machine learning is fraught with places where you're introducing biases that you didn't know about.” – Joshua Bloom
JOSH BLOOM: I'm Josh Bloom. Hi. Welcome to UC Berkeley. Welcome to Doe Library. Welcome to the Berkeley Institute for Data Science. And welcome to day two of Astro Hack Week. I hope you had a wonderful day yesterday. Today, we're going to talk about machine learning. I wanted to try to put machine learning into the context of an inference space that's certainly not rigorous, but one that hopefully helps you level set a bit relative to the types of things that you were hearing yesterday, and to give you a sense of how I want you to start thinking about bringing machine learning into your daily work.
What I should say from the outset is: start thinking about machine learning as just another tool in your toolkit for doing inference. This is not like, oh, now that I know machine learning, I'm not going to do Bayesian statistics, or I'm going to stop learning physics because, boy, machine learning is so awesomely data-driven that I don't need anything anymore. There's a time and a place for applying different tools. And what I'll teach you today is some of the at least Pythonic ways of approaching some astronomy types of questions around inference using machine learning.
But let's unpack this a little bit. I've got three axes. I've drawn them as orthogonal from each other, but for those that are reading ahead, you probably see that they're probably not that orthogonal. One is the statistical space, which is going in the left and right direction, Bayesian on one side and frequentist on the other. While I didn't see all the lectures yesterday, my prior knowledge of the two speakers suggests that you probably got beaten with the Bayesian stick pretty hard. As it turns out, most of the types of machine learning tools that we're going to wind up using are much more on the frequentist side of things. And those that, again, are practitioners of these various statistical frameworks for thinking about data and inference will often refer to the Bayesian-frequentist divide as really the wave-particle duality, that it's really both. And depending upon the situation, it's more helpful to think about data and inference in one place or another along that axis.
In the other axis, up and down, that's much more where you bring physics to the table or not. One is what you might call theory or hypothesis driven, where you've started from first principles, and you've derived a model that should lead to some set of observations. And that's the sort of hat that you wind up putting on. The other hat, which is the opposite end of the spectrum, is one where you say, I don't know anything about the system that I'm studying that I'd like to make inferences on, so I'm just going to throw the data into some framework and hope that I get the kinds of answers out to the types of questions that — because you've taken this hack day seminar — you've asked them correctly, in a well-posed way, although oftentimes, you'll find that people will be challenged in that respect.
And then in and out of the page here is the computer science view of this world, which is: are you doing the computation on a small enough amount of data that it can fit into RAM? Because if you can do that, then you get access to lots of really interesting algorithms that don't need to be that well parallelized, versus ones that are potentially out of core, where you're taking chunks of the data, and you're doing inference on chunks of the data at a time. And over time, you're building up a better and better inference about the questions that you're asking.
So this is the space that you should be thinking about. Many of you have been doing inference already, or you probably can localize the types of problems and approaches that you've taken in this space. And again, this is not going to necessarily be orthogonal from each other. I just wanted you to see the stats view, the computer science view, and the physics view of this. Yes?
STUDENT: Almost all astronomers have generally worked out of the board in the lower left, and Facebook is into the board on the upper right. But it's not clear to me that many of the occupants of this are actually filled.
JOSH BLOOM: No, that's true. It's like the fundamental plane of galaxies or the HR diagram. Just because you can draw temperature and luminosity doesn't mean that they're going to wind up populating everything equally. So there's definitely places of higher density versus lower density. That's a great exercise for us or for you to come up with some interesting problems that live in those different quadrants. I have got some thoughts on that, but I don't want to spend too much time. But you're absolutely right.
And that, in some sense, emphasizes one of the critical points, I think, of this whole morning, which is that you're going to choose your inference framework depending upon the types of questions you're asking and the type of data that you have access to. If you're at the Facebook level, you have access to, of order, several billion images. And so if you're trying to do insight into images and you can actually afford to get lots of labels, that is, answers for your previous images or videos, then you can bring a framework together that learns a lot from data, where you don't have to know much about, let's say, object detection. You can just throw it through some whole system that figures that out, like a deep learning network. Whereas in astronomy, oftentimes, if you're working right at the edge of signal to noise, you'll wind up not using machine learning because there, you would really benefit, essentially, pulling the signal out of the noise by bringing all your prior beliefs, both from a statistical perspective and also from a physical perspective.
Just a quick overview of what I'm hoping to cover today, and there may be a little bit more we can get into as well. What is machine learning? Obviously, we'll talk about that in all its gory details and cut as quickly as we can to the chase to try to do this in the most Pythonic way that we can. Go into two types of approaches and questions you'd ask of the data. One is regression, getting numerical values out of your data, essentially from a prediction perspective. We'll do a little bit of a runthrough on an existing notebook that I've created. And by the way, that notebook is up on the GitHub for Astro Hack Week, so if you want to follow along, feel free. And then classification, which is saying now I'm not trying to get some inference at the numerical level; I'm trying to understand which class an object winds up belonging to. And then we'll talk about how you actually improve your model. The first parts will be just introducing the topic, giving a little bit of theoretical motivation for those. But mostly, I want to get into what one needs to do as a practitioner. And then we'll start thinking about how to get these into production environments.
Most of what I'll be covering is what's called supervised learning. This is where you have a set of labels or answers on an existing set of data. And now you have a new corpus of data which is presumably taken and obtained and reduced in the same way as your original data, and then you want to ask the same sort of questions on that. So if you had a galaxy sample of, let's say, 1,000 galaxies that all had spectra, and you had your own classification scheme for what types of galaxies those were, and now you want to apply that to all the spectra from Sloan, this would be an example of a classification problem, but also a supervised one. So I'll make the distinction between unsupervised and supervised as I go along.
All right. So what is machine learning? The short answer is it's an offspring of statistics and computer science. The long answer is that it's a set of models which aim to learn something about a data set and apply that knowledge to new data. There's lots of different ways to unpack that, and there's lots of different ways to try to localize machine learning as something that has some boundaries to it within the broader context of artificial intelligence. What's started happening over the last couple of years is that there's been this almost Cambrian explosion of different algorithms and approaches, and commensurate with that are lots of books about each one of those and lots of tweets and people writing blogs about how awesome their algorithm is relative to other people's algorithms. One of the main difficulties, I think, not being a practitioner in machine learning and not having come through the ranks of having learned it from first principles — which I put most of us in that category — is that it can often be very noisy and very difficult to figure out what's hype and what's actually very useful.
The short answer on that point is don't believe the hype. But you probably should believe the hype at some level, because there's a lot of really powerful tools out there. And a big job that you're going to wind up having is trying to cut through that and trying to understand which are the tools that are useful and practical, and which are the ones that are highly specialized and really only need to be used in very specific cases. So it's the difference between a generic hammer that you would buy at Home Depot and some very specialized tweezer that can only open up a Nexus 6 or something. Very, very different things. One of them you absolutely need to have when you need it. The other one, though, is much more practicable and useful in your life.
Let me just step through these various components of what machine learning is, and then I'll break that out at an even higher level for you into different types. So it's using labels from training data to classify new objects. I've just given you a light curve. Is this a supernova or a nova? Domain experts can look at those and go, I know the answer. But if you want to do this at scale without any people in the loop, then that might be a good thing you'd want to use. I have a bunch of images. What's the galaxy type of these images, et cetera.
Learning the relationship between explanatory features and response variables to predict new data. What does that mean? There's a couple of terms in here that we need to get straight. Features are the variables. They're the x variable from what we saw yesterday. This is the input data. Oftentimes, we think of x in an abstract way. X could be an image. X could be a set of images taken over time. X can be metadata that's derived from those images. In the context of astronomy, it would be, let's say, photometry on a single object as a function of time. It can be photometry as a function of time on a given object, but it can also be the weather conditions at all the telescopes where that happened. This is all the data that could potentially get brought to bear on predicting this response variable, or what's sometimes called labels. And that's the y vector from yesterday.
Another utility of machine learning is its ability to discover natural clustering in data. Obviously, if you have lots of two-dimensional data and you plot it all up, your eye is going to wind up picking out clusters. You could draw circles around those and write papers about the objects in that circle. And that may actually get accepted to journals. And that's completely appropriate. But there's a more robust way to come up with the notion of clustering, even a two-dimensional space. There's obviously some parameters that define what it means to be a cluster or not that you'd have to decide based on your problem. But where machine learning tools wind up becoming really useful in this context is its ability to find clusters in very large dimensional spaces. If you have not two dimensions or three or eight, but 10,000 dimensions that you're looking at, and you're trying to find clusters that actually have meaning, you can use a bunch of different tools for that.
And similar to that idea is the idea of taking a very large dimensional set of data and reducing its dimensions down to just a few that are actually informative. A classic example of that in astronomy would be taking a bunch of galaxy spectra. While every single flux value at every single wavelength has its own meaning, and you can interpret that in a physical way, there's only a small set of representative spectra of galaxies. And every galaxy, in some sense, can be represented by that. If you've done principal component analysis on spectra, you've seen this, where if you take lots and lots of Sloan spectra at the same redshift, you can wind up picking out just a few representative samples. And then you can get an admixture, even in a linear sense, to build up most of the signal that you see in new galaxies. So that would be an example of taking a very large dimensional space — in this case, flux as a function of wavelength, potentially thousands of measurements — and condensing it down to just essentially a few weights on eigenspectra.
I should say, that particular example of PCA is where you start blurring the line of what is a machine learning technique and what's not — PCA is in the gray area there. Oftentimes, PCA — and this is a bit of an aside — will be used for doing dimensionality reduction on very large data sets with a large number of dimensions. And then the output of that dimensionality reduction would then be thrown into a machine learning framework. So in the context of, let's say, galaxy spectra, you might take the weights on the top 50 principal components to reproduce the spectrum that you saw and throw that into a more traditional machine learning classifier.
And last, because you can do clustering and you can find things that are like other things, the flip side of that is also true: you can find things that you haven't seen before. And anomaly detection is a very important thing in astronomy.
What I was going to try to unpack a little bit here — these are all the different components of what makes up machine learning. There's another way that I'd like you to also start thinking about what machine learning is and what it isn't, and that's in the types of questions that you wind up asking on the data. I'll call one type type 1, for lack of a better term, and that might be called exploratory. And the other one, which I'll call type 2, might be predictive and prescriptive. So why am I making that distinction? Exploratory in the context of, let's say, anomaly detection — that's a great example of what I'd call a type 1 activity with machine learning. You have a large amount of data. You want to find something you've never seen before. And then you want to look at that and say, wow, that's interesting, I should get more spectra of that object, let's keep going. And eventually, that could wind up leading to a paper.
It could also be that it's exploratory just visually. So this is a way of, let's say, taking very large numbers of dimensions and bringing it down to a smaller number of dimensions so that you could even visualize and play with, to get a better sense of the data and the low-dimensional structure that is effectively embedded in larger dimensions. I don't think of exploratory, in the context of lots of types of machine learning that's out there in the wild, as all that useful, because it doesn't oftentimes give you a tremendous amount of insight about what you need to do next. In the context of Facebook, the reason why they're building new frameworks for doing inference on videos and images is not because they want to learn about when somebody is happy, are they wearing this color or that color, and they're going to write a paper which is going to get accepted in some journal. They're doing it because they want to make predictions about what they can actually get out of that image. And then ultimately, they want to monetize the results of that.
And it's prescriptive in the sense that the best types, I'd say the more leading types of machine learning, are ones that say, not only do I understand what happened here in this data that you've just acquired, I'm now going to make suggestions about what new data needs to get obtained in the future to get to some objective function. And so while all of machine learning is an optimization — and as long as you can start couching your frameworks around optimization, at some level, you're doing some form of machine learning when you include data — I think the type 2 are very, very interesting.
In the company that I started, it's called Wise.io, we're working on more of the prescriptive kind, where we're looking at interactions between people, and we're making predictions about what they're saying and then making suggestions about what could happen next. What's nice about that type is that it winds up creating nice feedback loops, where if somebody takes an action that you didn't suggest, you wind up getting counterfactuals to your models. And then the models themselves can wind up improving over time.
I actually think while this is very exciting and there's lots of classes of machine learning frameworks that are being developed here, a lot of what we do in astronomy is much more on the exploratory side. Going back to this example of classifying large numbers of galaxies, that is saying, I'm going to build up a catalog of galaxies of this type. I've got a couple of exemplars from that space. I want to find more of those things. So you go and you find more of those things. And then it's really up to you to decide what to do with that. Not baked into the optimization technique traditionally is what you need to do next. So you're using it as a starting point to say, aha, I just found some interesting data that I didn't know existed there. And now I'm going to use that, because I understand the physics of what I'm looking at, to go off and do something interesting with it.
Is that clear on the distinction? By the way, it's obviously not a very hard and fast, black and white set of rules. I just wanted you to start seeing the two differences there. Any other questions about overall statements that I've made about what machine learning is? Yeah.
STUDENT: Maybe you said this and I missed it, but couldn't we say that detecting low-dimensional structure or lower-dimensional structure in high-dimensional data is all of these things?
JOSH BLOOM: The question was, couldn't you say that detecting low-dimensional structure in high-dimensional data is all of these things? I don't quite understand the last piece of that.
STUDENT: Well, all of the other stuff up there you could rephrase in a way that makes it look like detecting low or lower-dimensional structure of high-dimensional data.
JOSH BLOOM: Yeah, it's a good point. At some level, you've got a function on your input data, and you want to predict something else. This is a generic statement of making models on data to predict outputs. It could be that I take a very large dimensional data set, and my output is three numbers. And that's my way of reducing it down to something that I could potentially visualize and play with. It's not the same as saying, what class do I belong to? Just because I can take a 100-dimensional data set and then show that it's got some actual structure in a two-dimensional projection doesn't mean that I know what class of galaxy that is along those two dimensions. So you're right in that machine learning is a summarization of the data, and it is a compression, in some sense, of the data to get out an answer that you care about. But they're not identical. Any other questions? OK.
We're going to now start looking at the landscape of machine learning capabilities within Python. For those that haven't already, you can do a conda update scikit-learn, as you see up there. If you're using a different package manager distribution system, you can do a pip install scikit-learn. I believe the most current version is 17.1-1. It's probably fine if you have 17 and beyond. Some of the things in the notebook might wind up breaking if you're at 15 or 16. So I'll give you a second to do that.
While you're getting your environment up to speed here, what I will say is that there's a large number of packages that exist, even in the Python sphere, that do some parts of machine learning. There's mlpy. Orange, which has been around for quite a long time, is much more visual. Keras is a wrapper around a couple other machine learning frameworks, and particularly deep learning frameworks, and that gives you a very high-level view of deep learning. Nolearn is a competitor, at some level, with Keras. And then there's some of the lower level things, like TensorFlow, et cetera, that give you a little bit more direct access to some deep learning frameworks. And then astroML. There's a whole bunch of these things that are out there.
By far, I think it's fair to say scikit-learn is the de facto starting point for doing any machine learning in Python. And that's because of a number of things. One, it's very well-maintained. It's got something like 20,000 commits at this point, hundreds of committers. And they have a very well-documented and really sharp-elbowed API system. So you can wind up writing an entire pipeline against one type of machine learning model framework and then essentially rip and replace that out and stick a new one in once you wind up realizing you want to try something else, and not break a whole bunch of the rest of the code. Any questions about Python, machine learning frameworks in Python? Anyone else beg to differ? OK. This is a fairly complicated slide from the scikit-learn folks of what they call the cheat sheet for doing machine learning. I'm asking you questions about what type of data you have and what kind of approaches you're going to take, and then what's available within the scikit-learn world. The first one is, do you have more than 50 samples? If not, go get more data. So stop. That's not completely fair, because of course you can do inference on just 50 data points or 50 instances. Are you predicting a category? And do you know what the label is? If you're predicting a category, you're going to do some sort of classification. If you're predicting a response variable quantity, you're probably going to wind up doing some type of regression. Do you have large-dimensional data? Et cetera, et cetera.
I'm not going to flow through all of this, but what I will say is that once you've arrived in each one of these different boxes, what scikit-learn gives for you is a whole bunch of choices of the different algorithms that you can use. And one of my fairly strong statements that I'd be happy to be challenged on is that there's only a few algorithms out there now that are so battle-tested and so useful, that you probably don't need to learn nor try all these different ones once you've arrived in one of these different boxes here. That is to say, if you've got low-dimensional data and not a lot of it, some form of a linear model is probably fine and appropriate for your data. If you've got a medium amount of data, let's say at the gigabyte level or even the terabyte level, something like a group of decision trees, something called random forest, will almost always get you to the best answer possible. And then if you've got much more data than that, typically you'll wind up using some type of deep learning to get access to the best results.
It should also be clear, because all of you do do computation in one way or another, that the best answer, the one that essentially maximizes the objective function — that is, one that gives you, let's say, the best accuracy in a galaxy-star-quasar separator — isn't always the best one to use. So anyone want to come up with some arguments why the best accuracy algorithm might not be the one you'd want to use? David?
STUDENT: It might be outrageously computationally expensive.
JOSH BLOOM: It might be outrageously computationally expensive. And Dave talked about yesterday the trade-off that you have to make as scientists between what's the right thing to do and what's the tractable thing to do, because you need to publish. If it takes the age of the universe to get 2% better accuracy on a star-galaxy separator, and you have to burn all of your campus's resources to do that, you're probably not going to do that. So computational efficiency, efficacy, is an absolutely important one. What else? I'm going to save these. So why not use the most accurate? We said computational efficiency. Yeah.
STUDENT: You overfit to your input data.
JOSH BLOOM: Good — overfitting. Given a sufficiently large number of tuning parameters, I can fit a sufficiently arbitrarily large amount of data. But then that doesn't mean it has predictive power on the data that I haven't seen yet. So that's a great thing of what we want to protect against. Yeah.
STUDENT: Just building on that, accuracy might not be your goal.
JOSH BLOOM: Yep. So accuracy might not be the right metric. Indeed, your goal could be to minimize the number of false positives. A great example of that: if you've got a drone that's flying around, shooting down other drones, or airplanes shooting down other planes autonomously, you probably don't want to have any false positives, and you start shooting down your own people's drones. You want to shoot down somebody else's. I don't know, I'm coming up with examples. That's an example where accuracy is really important, but it's not the most important thing. There, you'd want to create a model that drives false positives down to zero. In the context of false negatives, if you're building a classifier on images as it streams off of telescopes, and you want to find things in the sky that are new and novel, you want to minimize your false negatives. You don't want to miss anything that's new and novel. But in that case, you also want to not find and say everything is interesting. So the two extremes of that, in the context of, let's say, a transient classifier that's happening in real time: one would be everything in the sky is interesting right now, and the other side of that would be nothing is interesting in the sky right now. And both of those are horrible systems in production.
What other reasons might you not want to use the most accurate? Is the most accurate the most informative? That is, if I have a classifier which just beats all the other classifiers, and then I look inside of it, and I say, I have no idea how it's combining the data — and yeah, I can see the math, but I don't get any insights out of that — for astronomers, that could be a pretty big deal. Unless you're working in a regime where the proof is in the pudding — that I don't care how you got the answer, just show me that, in a calibrated sense, you get essentially the best, most accurate curve or the perfect metric — very often, you're never going to be in that space, because you have to write a paper eventually about the things that you've built. And somebody's going to ask, well, why did you get to where you got? So they may not be all that informative.
All of these are reasons that you'd want to take a bit of a pause when you meet somebody on an airplane and they say, I work at startup x or a large corporation in the valley y, and we have the most accurate classifier on something. And you'd say, well, that's great, but it's like when you talk to an astronomer in a field that you don't know, and you say, what about magnetic fields? And they go, aha, magnetic fields. These are the things that you want to keep in mind when you're actually working with real frameworks, recognizing that overfitting isn't just what could potentially happen when you work with one model. You could wind up doing some sort of meta overfitting, where you choose a model that happens to look best across multiple different types of models. And that's because you've gone through effectively a multi-trials problem, and you've essentially glommed on to something that's best, but within the noise. Did you have a question?
STUDENT: Yes. Going back to that you wrote on the board that as predicted was your (INAUDIBLE). So out there in the tech industry, how much work goes into understanding the systems they work with, compared to just making models that can predict new data and enable them to place ads and (INAUDIBLE)?
JOSH BLOOM: So the question is, out there in the non-academic world in industry, for the type 2 kind of predictive sort, how much do they care how much they can unpack an understanding of what their models are producing? I think the answer varies. My guess is that in the ad world, in the ad placement world, it doesn't matter, because there again, the proof is in the pudding. If I get epsilon better than somebody else, that means millions or billions of dollars in accuracy. And boy, I can throw the whole kitchen sink at that because I can afford the computation to do it. Again, in the context of understanding whether you're happy or not, to suggest a different ad to you inside of Facebook, I suspect there are people within that organization that do want to understand what it is that their models are doing, and what type of data do they need to acquire, what type of data is informative, what's not, for the purposes of making a better prediction — but the end goal is not to actually understand why.
I think in the context of, let's say, financial models — let's say a risk assessment model — there are plenty of reasons, the most important one being regulatory, that when you deny a loan to somebody algorithmically, you have to say why you denied the loan. Just like with your FICO score, which is a terrible model of your ability to repay a $10 loan from your friend in two days, you can actually figure out exactly why you got the score you got. There's a regulatory requirement around that. There are regulatory requirements around making some types of inference models essentially completely transparent of how they got to that answer. I believe the EU just started passing a law, or they've already passed a law, that's starting to make some of those actually required to be transparent and then actually known by the public, of why they got the answers that they got.
I think the informative component of machine learning models is vastly understudied in the academic circles. Most people are focused entirely on, is it more accurate than somebody else's on the same data? Does it scale better by a couple of different metrics and something else? And very little of it has to do with interpretability. And I think that's a fantastically interesting open subject.
OK, so supervised learning. We're going to use a set of training pairs — so this is an x vector and a y outcome — to predict new y outcomes when we see new x's. In the context of regression, that's predicting a continuous variable, which will be our y, from an input set of features. Lots and lots of things and lots of different approaches exist within Python and scikit-learn, from linear regression to very nonlinear model capability.
All right. Why don't we jump into the notebook now. I'll mirror displays and see if we can get this to work. If you go into Astro Hack Week, Astro Hack Week 2016, and then you click on day two-machine-learning, you should see the latest version of the notebook. I'm going to do this and make this a little bit bigger. And we're going to do some astronomy examples. I adapted a previous notebook from Python 2 to Python 3, but I try to make it backward compatible.
Just show of hands, who's using Python 3? OK, so most of the people in the back. And then who's using Python 2? Most of the people in the front. That's interesting. I did some clustering in my head — which could be wrong, might not be right. That's actually, by the way, one of the problems that you generally will wind up having with unsupervised problems, is that how good you did tends to be very subjective. And with supervised problems, you can actually define a metric that you can test yourself against.
STUDENT: (INAUDIBLE) that result. That result was pretty good. You should take a photo of that and see if it's (INAUDIBLE).
JOSH BLOOM: OK. And now you can't change your result. And you can raise your hand twice, because you might be using both kernels.
STUDENT: Yeah, so who is a Python 2 user? OK, who is a Python 3 user?
JOSH BLOOM: Oh, it's more right than it is–
STUDENT: (INAUDIBLE) OK, great.
JOSH BLOOM: OK, good.
STUDENT: That's ridiculous. [LAUGHTER]
JOSH BLOOM: But how good is my pattern matching? Boy.
STUDENT: Going to find (INAUDIBLE).
JOSH BLOOM: Yeah. All right. We're going to work some examples using the Sloan Digital Sky Survey data. Again, this should work in both Python 2 and 3. I wanted to grab photometry that's corrected for extinction for 1,000 quasars that have known redshifts. As you can all probably do SQL in your head, you can see where this is all going and what this gets. But we're going to pull over 1,000 quasars from Sloan. And you should be able to do that locally if you just run and execute this — oh shit. I'm running it on demohub just because I don't want to use resources on my laptop. If you run this, you should get the data.
What does this data actually look like? Probably appropriate to look at it. It's not getting bigger. So here, we have the object ID, RA, dec, dereddened u-band photometry, dereddened g-band, et cetera, et cetera. Somewhere, we've got a — we don't have a redshift here, or we do have a redshift? Thought I pulled over redshift. One of these is redshift. We've got RA, dec, et cetera. OK, so that's the data. We just now pull over some other stuff you need — pandas, seaborn, matplotlib stuff to actually do plotting — and then pull over scikit-learn. I don't think I need to run this, but I'll just try it.
Does everyone know what demohub.jupyter.org is? It's a site that is built and maintained by the Jupyter collaboration, where one can have not ephemeral notebooks, which you can get from try.jupyter.org, but actually dedicated computation space. So you can do a git clone of a repository there and then pull in all of your stuff and then just have it there. And those computers, in principle, will run forever, and your servers will run forever. That's obviously not true, because it's maintained by a third party. I think it lives on top of Rackspace. But the whole machinery for building these things is actually now out and fairly robust, obviously open source. So for those of you that are teaching, it's a nice thing to be able to build a space for you and your students to all work together that's in the cloud. This notebook, obviously, I could be running from my own machine if I'd like, but I'm just running it from the cloud now. And that'll probably come back to haunt me because I just extolled its virtues.
So what are we going to do? We're going to import numpy, matplotlib, pandas, seaborn, et cetera. And we'll take a look at that data. Just doing a head on the first 1,000 quasars. And you see all the stuff that I got out of that. Oh, here's the z. Spectroscopic redshift is spec_z. And here's the class. Now, we're not going to just do machine learning on it. One of the things that you have to know is that you have to look at your data before you start throwing it into these frameworks. Just because these frameworks are nice and easy to use doesn't mean that you're allowed to get away with not actually looking at the data. So let's look at the data.
What are we going to do? We're going to pull only some of those columns out, because we're going to wind up doing a regression problem or, later on, a classification problem. I'm going to drop something like right ascension and declination, because I want to not care where you are in the sky if I want to decide if you're a quasar or a galaxy or a star, for instance. This is a great example of bringing prior knowledge to this problem. It should be the case that if you build a nice classifier or regressor that included right ascension and declination in the data set, that you would find that those features, those values of x in that vector, have no informative value and no meaning. But because they're in there, and because they're going to wind up not being distinguished between something else that may be much more informative, putting on our domain-driven hat, we have to allow ourselves to use domain knowledge just to make our models more understandable, more tractable, and actually have some meaning.
So we're going to pull out the object ID. We're going to pull out the spectroscopic redshift. We're going to pull out the different colors, et cetera. And then we're going to make features, which is just a copy of that whole data frame. And then we're going to get redshift. We're going to try to build a predictor of redshifts using just photometry data. And we're going to delete, out of our feature set, the answer, because that would be pretty bad if we had the label in the thing that we're trying to use to then predict our y values. And then let's see what the result is of that.
So this is our features. And now we've got a smaller number: one, two, three, four, five, six, seven, eight, nine. And some of these features, by the way, are the difference between aperture photometry and Petrosian photometry, which could be informative for the size and the shape of the galaxy at some level. I'm not taking all the other metrics that are obviously also available. I'm just trying to do something simple here. So we've reduced this down to a nine-dimensional input y vector. And for those that are used to pandas, I've created an index around the object ID.
And now, let's plot the histogram of the output variables. And what you see right away is that there's very little numbers of quasars beyond redshift of 2 in the Sloan catalog, or at least of the ones that I wind up picking out. So there's a pretty big cut around 2. So you can imagine that when I create a regressor, if I don't give it any other information, it's going to try to do a very good job around the places where there's lots of data. And it might not know that actually, what I really care about is finding large redshift quasars. So again, this would be how you would wind up constructing the model. You'd have to think about that.
So let's plot the data against itself, pairwise, every single feature against every other feature. Let's just see what looks informative. And we'll colorize that by redshift. This takes a little while to produce here, so give it a second. Sorry to Stefan, who's sitting over there, because I'm using Jet and not your color map. There we go. OK. So there it is. And it's going to be very hard to see here. We've got dereddened r, and it's just a histogram of that down the diagonal. And you can see indeed that there is correlation, heavy correlation, between the difference in the u-band between the two different photometry methods and the g-bands and the i-bands and the z-bands. Those ought to be correlated with each other, and indeed, they are.
You also see something kind of crazy here. What is the dereddened r-band magnitude? It looks like it's mostly around zero, but some of them are around minus 10,000. This is like somebody setting off a supernova in your eye or something. So that, obviously, isn't valid data. And we just realize that there's something crazy with the dereddened r-band magnitude. Some of these actually are not detected in r-band, and we can decide what we want to do with those. Again, here, we're not even doing any machine learning. We're doing what's called feature engineering. We're going through a process of understanding our data. We've articulated the question that we're going to wind up asking: can we predict redshift from just photometry on quasars?
And what we'll now do is we'll say, look, I'm just not going to care about all those ones that don't have r-band magnitudes. Maybe those ones would be really easy to find what their redshifts are. But just for the sake of this demo, if the color is sensible and it's actually a sensible magnitude — that is, it's not this minus 999 — we're going to keep it. Otherwise, we're going to get rid of it. So now we can look at that matrix again and make sure that we get what looks like pretty sensible answers. And then once we're done with that, we'll wind up saving that data to a CSV file so we can use it later on and bring it back in. This also takes a little while. I guess these computational clusters we're running on are not all that beefy. Maybe I should have used my laptop, but we'll give it another second here. There we go.
So it may be very hard to see from the back of the room, but you can start to see, in some of this data here — let me see if I can zoom in even further. No, that doesn't help. You can start to see, in some of this data here, that there's different colors. There's some green and there's some blue. And it starts to look like there may be some separability in that. So this should give you a little bit of hope that you may actually wind up being able to figure out the different redshifts in this large dimensional space. But it's very clear that it's not like in two axes, we get perfect separability of high redshift versus low redshift, for instance. So this is actually a pretty hard problem to do. And those that have worked on things like this know that it's very hard to do it.
All right. But we'll say and declare, for the purposes of this demo, that we're done with the pre-processing steps. I'd still call this part of the machine learning pipeline that we'd have to build if we're going to do this for real. But now we're ready to do some basic model fitting. To do basic model fitting, we need to create a training set, and then we're going to create what's called a testing set. We're going to hold back some data on the side, and we're going to take our model and predict against data whose answers we already know. And then we're going to try to see how well we did. So we're going to create our x vector and our y vector. The x vector is going to be features, and y vector is our answers. We've got 9,988 examples of data that meets our filter criterion in nine dimensions, and we've got an output of the same size. So we're off to the races. Just for the purposes of this demo, we're going to choose half of the data for training and half of the data for testing. What did I just do, other than just create a train set and a test set? What did I just potentially introduce in terms of biases into the thing that we're going to wind up building? Yeah.
STUDENT: (INAUDIBLE)
JOSH BLOOM: That's right. So it's not a random selection of the data points. I assume that what came back from Sloan was not something ordered by magnitude or ordered by location in the sky, which could have some very subtle impact on whether you're in the supergalactic plane or not in the supergalactic plane. We just introduced a whole bunch of biases in the way I even created this train and test set. A better way to do it, or at least another way to do it, would be to create a random test set. We'll see that Python has a whole set of methods that allow us to, without having to do what we're doing right now, implicitly create train and test sets as you wind up building up the model. But good. So really important to bear in mind the assumptions that you're making about the data as you wind up building up these models.
So we can build a linear regressor on that. Typically, what you'll wind up doing is instantiating one of these different models, and you call it clf. And you can wind up doing a tab completion to see all the different things that can be done on this model fitter. The first one obviously is to run fit. Once it's been fit, we can get the parameters of the fit, et cetera, et cetera. And we can apply that, if we want to, to new data.
So let's fit the data. For all the learners, they will all have a dot fit method associated with them. And so it should be where you can swap in and out these different learners and still apply it to the same format of the data. So we'll train on the x data, given the y outputs. And we got a result here. We're going to wind up predicting on our test data, and then we're going to see how well we did. We're going to get the mean squared error relative to our test data. Visually, we can see what that actually looks like. And you see that this is pretty crappy. Yet at some level, it's actually pretty good, because it did the best job it could for a linear regressor on nine dimensions in fitting most of the data. But you see that it was a massive underfit of the high redshift results.
Here's what we would do if we just took the average of the training data. So mean squared error of 0.65. And if we'd just chosen the average, we probably would have gotten about the same result. So this is not a very good regressor. And you can see some things here from scikit-learn, where we just don't have to build up the notion of what a mean squared error is ourself. They have all these different scoring functions.
Before I jump into other classifiers, what I'll do is jump back into the lecture notes, and I'll introduce some of the other — I won't call it theoretical motivations, but at least visual motivations for these other types of classifiers. But anyway, I wanted you to just see at the very beginning here that we've got reasonable end-to-end of starting off with raw data, and then we're building a model, and then we're able to use that model on new data. It just turns out to be a really crappy model.
STUDENT: I have a quick question.
JOSH BLOOM: Yes, question.
STUDENT: When you do this mapping, what does it actually do? Can you tell us?
JOSH BLOOM: What does the fit actually do when we do the matching? Give me one second. Just put something up here. Sorry. The question was, what are we actually doing when we're doing the modeling?
STUDENT: The fitting — because you use the first data set, and you fit a regression model to this data set and then use that model to the–
JOSH BLOOM: So the question is what's happening under the hood. When you run dot fit, the machinery of scikit-learn is going through, looking at your input data, looking at your outputs, and coming up with a linear model that winds up predicting this. So it's coming up with — however the linear model is created and run within that; there's logistic regression, there's all these other ways of doing modeling — it's just coming up with essentially weights on the data. What are the nine numbers I need to multiply to get out another number? It's just coming up with that. And it's storing those results so that when you then apply that to new data, you wind up getting out essentially what the result should be. And because you know the results on the test data, you can compare them directly.
What I wanted to show you is this thing called k-nearest neighbors. K-nearest neighbors is represented right here, which is to say I have a new data point. This is showing you a data point in two-dimensional space, so I'm not showing you what the axes are here. But I have a new data point, which is green, and I want to know what its value is, its result is. In a classification sense, you might just say, well, if I have a parameter which is called a hyperparameter, k equals 3, I'll just take the 3 around me. And if I know the two classes are of class blue and class red, then I can wind up saying, well, probably with a 66% probability, this class of green dot should actually belong to the red triangle class. If k is equal to 5, you may wind up saying that it's equal to the blue class. So this is a tuning parameter of a model like this.
In the context of regression, what this is doing is saying, let me take the y value variable outputs of those three things inside of that first box. And I'm going to, let's say, average those together. Or I could take the median of those. And you can define exactly how it winds up doing the combination. And then that should be your answer. And this works in any large number of dimensions. But as you can imagine, you get into very weird places, where you're in 1,000-dimensional space, and really, nothing is near you. Almost by definition, nothing will be exactly near you unless you have a repeat of the instance from your training set. And so the notion of what that distance is becomes really strange.
One of the things that's really weird about k-nearest neighbors and linear regressors and lots of the types of machine learning models that are out there is that it assumes something implicit about the data. Does anyone want to guess what that is, based on some of the words I just used? I've got a nine-dimensional space in the problem I just used. What are the units of all those nine dimensions? Those happen to be magnitudes, right? If we go back to the notebook here, I've got magnitudes. But look at this. Most of these are color differences, so their values are pretty close to zero. But I've got one, which is the apparent magnitude of the object, and these numbers are nowhere near zero. When you build a linear regressor, there's no difference between this value and this value here. And when you're doing k-nearest neighbors and you're trying to find out the distance between you and the nearest object, the probability of having a distance of 0.4 in u minus g is much different than the probability of having 0.4 of the difference in the dereddened r-band magnitude. All this is saying is that many of these different classifiers and regressors that we wind up using have this implicit assumption that the data across all the dimensions are of equal importance and have an equal notion of distance, in a Euclidean sense. And so when we build a distance metric implicitly in some of these machine learning classifiers and regressors, we're implicitly assuming some metric in that space.
What would happen if I added another column here that was something that was potentially qualitative, like near a galaxy or not near a galaxy, or near the galactic plane or not near the galactic plane? How would a linear regressor deal with that? One possibility, of course, because now you're dealing with categorical features, is you could turn those into zeros and ones and say if you have a galaxy near you, we're going to call that zero, and if you don't, we're going to call that one. So then you have binary variables as some of your input, and then you have continuous variables as some of your input as well. Many of the original classifiers used in machine learning don't do all that well in that context.
Let me now do a k-nearest neighbors and see what happens, and see if we actually get an improvement in our model here. We already saw the result of linear regression, and we saw that the mean squared error was not much better at all than what we got from just essentially choosing the mean. So here, we're going to wind up now dealing with this metric space and doing something that tries to put all the different nine-dimensional parameters on the same footing. And that's where we're going to wind up scaling them. We're going to wind up scaling them all to approximately the same distribution, so they'll have, I believe, a mean of 0, and they'll have standard deviation of 1. And this scaler winds up remembering what it did to the input data before it wound up running a k-nearest neighbors regressor.
So I'm going to use k-nearest neighbors of 10. I'm going to take the nearest 10 objects in my nine-dimensional space, and then we're going to use whatever the defaults are from scikit-learn. And we're going to wind up fitting. You see that fit was pretty quick. It assumed a Minkowski metric, which is fine because I already did the pre-scaling. And we're going to see what our mean squared error is, which is much, much improved. So now our mean squared error is 0.23, as opposed to 0.6-something. And you see that I still have sort of the same structure I had before, so that I'm missing some of the high redshift stuff. So effectively, if you think about it, these high redshift things at redshift of 6 were assumed to be off by 4. It assumed all the high redshift stuff was at redshifts of around 2, which is near the mean. And you see, I got a result here which is much better than what I had before. But I have a hyperparameter. So is that result better because I chose k equals 10 or k equals 5? Let's see the results out. That wound up improving things a little bit. Let's try k equals 1, so we'll just take literally the nearest neighbor. You see the result is much worse. And you can start to see some structure in the data here, which is a little bit scary. I don't know, let's try 20. The results are all about consistent with each other.
So here's an example now where we had to do some pre-processing. If we didn't do the pre-processing — I think we can try this. What would happen if we didn't pre-process that? We already pre-processed x. I don't want to go back up and reprocess all that. But it'd be worth just trying it for yourself, what happens if you don't do this pre-processing step here. OK. So here's a pretty good example. I think I played around with this for a while, and I found that k-nearest neighbors of 5 gave me approximately the highest accuracy answer. Any questions about that, k-nearest neighbors? I'm going to jump back to another type of classifier and regressor called random forests in just a second. But if there are no questions, then I'll jump over.
By the way, when you split up the space with k equals 1 — this is, again, in two dimensions, so I've got three classes, red, green, and blue — this is the space that, if you now sample densely across your learned space, you get a very different classifier. That is, for the same point that shows up — just randomly choose a place in your head there — you often wind up getting a different answer depending upon what k value you wind up using. There's a nice thing on the web, by the way, I wanted to show you, where you can see — I don't know, these are just support vector — all right, I'll come back to this. It's support vector machines and random forests. I'll be able to show you visually and interact with that. Let me show you now a little bit on decision trees. Yes.
STUDENT: (INAUDIBLE) comment earlier, but could you actually talk about your favorite clustering algorithm (INAUDIBLE)?
JOSH BLOOM: What's my favorite clustering algorithm? I often don't do clustering. I feel like the most interesting man in the universe. I don't do clustering, but when I do, I drink Dos Equis, and I do k-nearest neighbors. But those tend not to do great in very large dimensions. So there's something called DBSCAN, which tends to be quite good. And oftentimes, when I'm doing clustering, it's because I really just want to see low-dimensional embeddings. So I'll often use some manifold tricks that I'll show you towards the end of the lecture to do that. And then most of the time, again, because I'm generally in the business of doing lots of unsupervised problems, I'll wind up trying to couch the problem in a supervised way. And then I won't use k-nearest neighbors. Questions?
STUDENT: What if you have, for example, periodic features? So right ascension of 0 being very close to right ascension of minus 1, which is the same as et cetera, et cetera, et cetera.
JOSH BLOOM: So the question's about periodic features — right ascension being very, very close because you loop over on the celestial sphere, for instance. That's a great example of where domain knowledge needs to come in. If you don't care about location on the sky, throw that out of your feature set. But if you do, because that's going to wind up being informative, put it in terms that you believe will be most useful for the output that you care about. So if you want to know, for instance, that this thing that you just found is an asteroid, rather than give it RA and dec, you probably want to give it ecliptic latitude and longitude, or ecliptic latitude only. So that would be an example of taking that data, changing it to get closer to what you think is going to wind up being informative. The probability of building a — I shouldn't use the word probability so easily — the chances of one building an amazing classifier which winds up learning how to transform RA and dec into galactic latitude is more or less zero. You're not going to build something that's going to wind up inferring the math behind that. If you already know the math, by all means, you should be entitled to actually apply that in what you would generally call the featurization step. As most people who do machine learning, both in industry and academia, know, when you're applying it to real data and asking real questions of it, you're spending almost all your time in featurization. That munging that we just did, which I breezed through really quickly and was just really simple filtering and stuff — that's the kind of thing you're going to wind up spending a whole bunch of time on, because it's the place you're not only introducing more informative features, it's also the place you're starting to protect yourself against bias. There was another question over there.
STUDENT: So you need to standardize the data in the second round (INAUDIBLE) and I was wondering how (INAUDIBLE) and how much that affects (INAUDIBLE).
JOSH BLOOM: So the question was, we just scaled the data as a pre-processing step before we did the k-nearest neighbors — how much do the parameters of that scaling wind up impacting the results? The short answer is, I don't know. We just have to try it. And if we were going to build a real paper around this or a real result around it, we would spend a lot of time trying to understand how sensitive our results are to the different scaling techniques. And it may be that we want to get our results very sensitive to the scaling techniques, because we want to highly optimize on the data we have. The flip side of that could be maybe I want to build a star-galaxy separator using k-nearest neighbors that don't just work on Sloan data, but work on other types of data, in which case we'd try other types of data after we finish that whole exercise and convince ourselves that we got the types of accuracies that we were expecting.
So we really don't know what the impacts are of the choices that we're making until we actually explore that and try that. And this is the crux of why you all would be better at building a star-galaxy-quasar separator than if I handed this to the most advanced data scientist who worked at pick your favorite company, Facebook. They might know a lot more machine learning than you. They might know a whole lot about all the intricacies of pre-processing of data. But because they don't understand the domain, they wouldn't do as well. And they'd potentially overfit without even knowing it. Any other questions?
All right. Let me jump now to decision trees and random forests. I'm going to present this section mostly in the context of classification, but it also has implications for regression as well. Before I do that, let me just show you a very classic data set that people wind up using in machine learning. This is called the iris data set. It comes pre-packaged with scikit-learn. And here, you've got essentially three dimensions — sorry, four dimensions — measurements on these various different flowers, and three different classes of what type of iris these flowers are. Has anyone not seen this before? OK, so I'll breeze through it because only a few of you haven't seen it. This is the thing that you will typically apply your classifier or regressor against to see how well it performs. Here, you can visually see that I can make cuts on the data that allow me to separate out different classes very easily.
For instance, in sepal length and petal width, if I made a cut in the top right quadrant over there, up and down vertically, around 0.5, I would be able to perfectly remove all the classes of red. And then I'd be left with needing to make some other cuts in another dimension that wind up giving me the separation between the green and the blue. So this is in a classification sense. If I now wanted to say, for a regression sense, given the values of, let's say, the first three parameters, what's the fourth parameter? Can you actually predict that? That's another question you might ask of the data.
But let me step back to decision trees and how it views that. This is a decision tree that's built on the iris data set. All it's doing when you build decision trees is you're building essentially these branching points, where you wind up, in the case of prediction, when you wind up having a new input y variable, at the very top level, we're saying, is the petal length greater than 2.4? And if it is, go to the right. If it's not, go to the left. And if it's in the left, you see what's called the terminal node up there. Do I have a pointer? Yeah, so this is what's called a terminal node right up here. But if not, I'm going to go down here, and I'm going to ask another question. Well, is the petal width greater than 1.65? And if it is, go down to the right. And if it's less than 5.5 in the petal length, go down here, et cetera, et cetera. And eventually, I get to another terminal node for the given instance of what I was just given that tells me that this was of type Iris virginica. So that's classification. In the context of regression, the idea here would be to go all the way down, and you would wind up just getting the value of one of the examples from the training set.
How you actually build up these different trees is kind of interesting, and that's where the real interesting math lies. But the way to think about it is you wind up trying to figure out how to create the most amount of separability between the classes that you know about, or if you're doing a regression problem, you try to get the most information out by making one of these splits. There are always splits at every single node that happen for one of the variables in your feature set. And how those variables are chosen as you progress down the node is where all the interesting math winds up coming in.
But you can see here that I can build a perfect classifier. In the x2 dimension, if it's greater than 1, then go down this path. If not, I've just gotten the green class. That's a perfect classifier on that. And then if x1 is greater than 1, then I will wind up splitting here, and I can classify between blue and red. How you actually do that is beyond the scope of this tutorial, but there's a couple of different ways in which one does it. Most of it — the way to think about it is you're trying to gain the most information by making a split. What you will typically do, when you're building up one of these trees, is you'll randomly choose a feature to try out of your large-dimensional feature set. And then you'll move effectively a pointer around for all the different values that you could split on that. And you'll stop the pointer where you wind up, on the other side of it, getting less information by cutting. And the information gain, or an entropy notion from a Shannon information theoretic perspective, are just simply defined these ways over different numbers of classes. Or you can define a similar type of thing when you're doing regression.
But depending upon exactly which feature gives you the best gain, you will just say, OK, for this node, I'm just going to cut on that one. And then as you go farther down into this tree, you wind up taking all the examples from the training set that made that cut, and you'll keep on going. And then you'll have a stopping criterion of when you decide, I'm actually done, and get to a terminal node. Oftentimes, you'll stop when you get to one single instance of your training data. But many times, what that means is that you've got lots and lots and lots of cuts to get all the way down to isolate every single training example. So one of the parameters of these decision trees will be, what's the minimum number of objects that I need to have in my node before I wind up stopping, or the maximum number before I stop? You can also imagine, as you wind up building up these trees, that you would actually want to stop when you get to a fixed length.
Now, it turns out that trees are really good for classification and regression. And what I'll show you here is that given a classifier, I can make a decision, when I put down some new data, exactly where I'm going to wind up showing up in this space. But what's interesting is that if I show up in this space here for, let's say, new data — that is, given a new x2 and x1 — not only do I get a good classifier (so if I'm in this space, the answer is green), I also get some notion of probability of which class you belong to. So here, I had one example of blue, and I had a couple of examples of red. If I'm somewhere in here — and that's the spaces I partitioned using my decision tree — then I wind up actually being able to get some notion of probability, which can turn out to be very, very helpful and useful.
Now I can show you a demo. Here is a demo of decision trees. This is essentially the result of having built not just one tree, but multiple trees. And when you build multiple trees, those are called large numbers of trees, or random forest. And it's got some other names associated with them. Now, if we've got a class of red in our training data and a class of green in our training data, this is the space as it's partitioned by not just one tree: as I put in essentially a new value here, which space do I belong to? You can also see that there's some notion of a probability. So out in here, I'm not sure whether I belong to green or red. And in fact, in this case, you're going to get an answer of what the probability is. But because you're so far away from your original data, this might be an example where you give pause and say, I've got a probability out, I got an answer, but I'm not sure if it's right, because I'm far away from my original data. That's another thing.
If I now plop down another data point — so I plop down red and build a new set of trees, and now took an average over a bunch of these different trees to get what the result is — you see that I wind up changing the space. If I put one right inside of here, I wind up changing the space; put here, et cetera. So if I want to do a green data point, shift-mouse click. You notice here, I'm not changing anything except I'm making this quadrant very, very strong for green. If I put a green right in here, things get pretty ugly. So I just drew a little red box around this one red point.
And you see that it's not just a fixed number of parameters. If I change the number of parameters of my model, I wind up getting a different space out of that. So if I have a few number of trees — in this case, I've got one tree — this is what one tree would wind up splitting up the space to be. If I have lots of trees, I split up the space — so this is 189 trees — I split up the space like this. The number of trees you have may actually wind up being part of your protection that you have against overfitting. If you have lots and lots of trees, in principle, you could wind up overfitting on the data because you've just seen it too many times. Random forest does a pretty good job of protecting yourself against that, not introducing too many biases with it. I can change the depth of how deep the trees can go, et cetera, et cetera. Any questions about what I'm showing you here? Yeah.
STUDENT: Is there a notion of error in the parameter? So the red and green points are next to each other. What if one of those has an error that (INAUDIBLE)?
JOSH BLOOM: Good point. So is there a notion of errors in the data? Yes, because there are astronomers in the room, and we know that our measurements aren't perfect. And there's noise in the data, and we have biases and errors in our measurement. But no, in the context that random forest doesn't take those into account explicitly. One could take it into account by making the error on, let's say, the r-band magnitude another feature. But what you're not doing is telling random forest — or any other of these classifiers, by the way — feature two is the error on feature one. It would just potentially wind up learning that if there's a large number in feature two, then I should discount the values in feature one. One of my big knocks against a lot of the academic work in machine learning is that very little of it contemplates uncertainty in the input data. Likewise, very little of it contemplates the fact that you might have mislabeled data. So in this example here, our original data-taking may tell us that this value is red, but the class is actually green. And so we just messed up our classifier locally around that point. So in a Bayesian sense of having parameters of some probability that every data point is incorrect, you might want to actually take that into account in your classifier — we don't do that. The only way to really do it is with some notion of a brute force, where you would simulate what your data sets would look like under different types of assumptions, so you could jiggle around each one of these different data points. It's mildly appropriate to create lots and lots of instantiations of your data by randomly sampling what you observed from what you believe to be your noise properties of what you observed, and then using that to build up what is effectively a much larger data set, to build up a potentially more robust classifier. But that's a bit of a cop out. Again, most of the machine learning that you're seeing here and as done in the real world is very frequentist. And we don't typically have a good way of incorporating errors from the data into it. Any other questions?
STUDENT: (INAUDIBLE)
JOSH BLOOM: So the question is, given that we know that we're changing the whole model space as I change the hyperparameters of the model, how do we go about choosing the right hyperparameters, the most appropriate ones? That's what's called model selection, a hyperparameter optimization. That's what we're going to do a little bit later on. It's somewhat of a black art. There's machinery within scikit-learn to help you do that and run all of it. But in the back of your head, you should always be thinking, as you're doing these model selections and hyperparameter optimizations, am I just overfitting on all of this because I just happen to choose the right one that gives me the right answer or the best, most accurate answer? But we'll talk a little bit more about that. It's a good thing that you're thinking about it.
I want to be cognizant of the time. When should we be giving some thought to stop?
STUDENT: About halfway.
JOSH BLOOM: All right. Why don't I go through, and I'll just finish up on regression, and then we'll come back and do a classification example. And then we'll do a little bit of anomaly detection after that. All right. Let's do some random forest. You instantiate a random forest regressor. Here, we're going to build not one tree, but 100 trees. We're going to create a criterion by which we decide on the different splits at every single one of the nodes — this is just going to be mean squared error. Minimum number of samples in a leaf will be one. These are different parameters that you're going to have to get to know if you're using random forest. So I'm running that now. You see it takes a little bit longer than the other ones did. And it just finished. So now, we can plot up the results. And you see that my mean squared error is now lower than all the other ones that I had before.
One of the major problems of random forest and these nonlinear learners — because it should be obvious that this is a highly nonlinear learner. I'm not partitioning the space with just a bunch of different hyperplanes. I'm actually making very, very nonlinear cuts when it comes to variable interdependency. And it's the job of a good learner to figure out the relationships implicitly between these various different features. But what should also be clear about random forest is that it allows you to make cuts on the data where you don't have to do this pre-processing step. And for that reason alone, it's why I tend to stay away from the ones where you do the linear regressors or the k-nearest neighbors, where you have to effectively define the metric space ad hoc.
In this case, because every time you wind up making a cut on the data, you're deciding to go left or right, you're deciding it based on the units of just that column of that feature. So if I'm cutting on dereddened r-band magnitude, I'm making a decision based on, what's your dereddened r-band magnitude at that split point? I have no reference at that split point to what the other values are in the other features. Whereas when I make a decision around how many neighbors do I go out to, you have to make a decision in a large-dimensional space, which means that all the features have to be effectively scaled correctly with respect to each other. One of the beautiful things about a random forest is that it doesn't require you to have to do that. And so it takes away one of those really ugly places where one, as a data scientist or a physical scientist, would have to imbue a lot of knowledge to get it right with some meaning. So now we have some questions about model selection — choosing the best model. We have this notion of what's called cross-validation. What I just did is I built a very poor man's cross-validator, where I took half the data, and I essentially took it out of my model building. And I built a model, and then I applied it to the other half just to see how well I did. But what you can imagine one can do is instead build a model on, let's say, 80% of the data, and then predict on the other 20%. And then move your window to another 20% of the data, and build a model on the first 20% and the last 60%, and apply it to your second batch of 20%. And now what you have is a prediction on all of your data, where you've now not used half the data, but now you've used 80% of the data. So in principle, your model should be better. And as you wind up tuning your model, you can imagine making those windows smaller and smaller. So this is called cross-validation.
Just visually, I'll show you what that looks like. This is what's called five-fold cross-validation, as depicted here. This is just a two-class problem, but it works as well on regression problems. You have a test set and a train set. And so you take the first 20%. You hold that out. You build a classifier on this portion here, this 80%. And then you apply it here. And then even though you're applying it on a small amount of data, you store the results of that application. And then when you're done, you basically have built up a whole bunch of predictions on data that weren't part of the training set. And as long as these models are trained on the same hyperparameters, as you go to each of these different so-called folds, you wind up getting a very accurate measure of what it would be like to have a model that's built on all of the data and applied to another set, which is out here.
By the way, oftentimes what people will do is not just create a test set and a train set. They'll also build what's called a validation set. In that case, let's say you have 120 examples. You take 20%, and you put it fully aside. And you say, I'm not going to do anything with that, other than to say, when I'm ready to write the paper, I'm going to apply my best final model to that and say, that's what I expect my results to be. And then you can do these clever examples with train and test to help you get a better and better model. But it's considered very good practice to hold out even more data in this validation set and not use it at all to help you train your model, because as soon as you start building models and changing hyperparameters based on the results of what you got on here, even though these data are, in principle, completely held out from here, there is some implicit crosstalk between the train data and the test data. So there's lots of different ways to do these different cross-validation holdouts, as it's called. The one that's often useful in the context of classification is something called stratified k-fold, where you're saying, for a large number of classes — let's say I've got 12 different Hubble types of galaxies, and I want to classify something in — if I have a class of galaxies that is a minority number, let's say there's only 10 S0's in my galaxy classifier in my training set, but there's 1,000 of some other type, then you can imagine that a classifier that's built correctly on the data, a little bit like we saw before with the redshift regressor, would try to fit most of the data. And so there's a nice thing about this k-fold, where it winds up holding out the appropriate amount of data for the given distribution of the classes in your training set. Don't worry if that's a little bit beyond what you've thought about for now. I just wanted you to hear it and think about it a little bit for now.
All right. So we're just going to build a little helper function, which is going to tell us what our cross-validation score is. We're going to do our cross-validation using k-fold cross-validation. We're going to do essentially 20% holdout, 80% prediction. Oh, wait. What? Oh, I forgot to run this. I get an answer there. We'll do cross-validation k-fold with 10 using random shuffling. This gets to one of the questions earlier about the potential bias that we introduced by holding off the first half of the data versus the second half. This will actually shuffle it for us.
And then we can do it on other classifiers. So this is clf2 — this was the random forest classifier. You see this takes a lot longer. This is building up a more intricate model. Let's see how well that does. Any questions while this is finishing up? Is this running on a Raspberry Pi? Is Rackspace all Raspberry Pis? I guess not. Well, we'll see. So what is this doing? It's actually doing a three-fold cross-validation, which means it's building three random forest models. And you see that our mean score winds up getting better. This, by the way, is different than the mean squared error, so this number should be larger.
All right. I'm going to stop here. We'll take about a 10-minute break or so. We'll come back and do some classification. We'll do some hyperparameter optimization. We'll do some parallel computing. And then we'll try to do a little bit of anomaly detection as well. I'm happy to take questions during the break if people want to come up.
Let's start up again here. I've been making reference to classification problems. Before we jump into that, I will just formally describe what I mean. Predicting the discrete class — that is, are you in class A, B, C, or D through n — of an object from an input set of features. So instead of it being a continuous variable that we're going to wind up predicting, we're going to try to decide which class you want it belonging to. And of course, you could turn a regression problem into a classification problem, where if you're above the value of 1, then you belong to this class; if you're less than this, you belong to that class. And vice versa. You could take a classification data set and say, if you're in this class A, that means that, from a regression perspective, I'm going to put you somewhere in the bucket between 0 and 1. But the closer you can get to defining your problem that's consistent with the question that you're asking and the data that you have, obviously, the better. Anyway, if I've got a bunch of instances — these are examples of x vectors and pairs of y, which are the labels, say 150 of them — and I want to predict these three different classes, that's actually what defines the iris data set.
So lots of different ways to do classification within scikit-learn. You can turn a logistic regressor into a classifier. K-nearest neighbors classifier — we saw that a little bit already. LDA, Latent Dirichlet Allocation, Naive Bayes, support vector machines, classification trees, random forests, blah, blah, blah. So there's a whole kitchen sink, like there was in the regression side, for you to make use of. And then there's also something which is also quite useful called feature selection, which is an automated way of deciding which features are informative and then rebuilding models based on what's informative and what's not informative.
Because I haven't shown it before, why don't I just quickly introduce yet another way of doing regression and classification: support vector machines. Typically, you think about these in the context of 0, 1 classification. So do you belong to the red box class or the blue circle class? And support vector machines is a mathematical model for building up what are essentially called hyperplanes — because we're showing this in two dimensions, it looks like a line. But in multiple dimensions, you want to create these nice places to separate.
So in the iris example, there's obviously some nice places that you would wind up separating between the red and the green. The question is where, and how would you build those lines? Obviously, I could build a line that cuts right close to the red and say everything to the left of this line in data we haven't seen yet is going to be of red class, and everything to the right will be of the other classes. But is that the optimal place to put it? The whole idea of support vector machines is that you try to get these hyperplanes in between the data, where you create what's basically the largest space right here. So this is the largest support, as it were. So this A is the optimal line to separate this data set and this data set. Even though B also does a perfect separation, it gets very, very close to this value here and gets very, very close to that value here. So it has less support. The larger the support you can get, the better. This works really well, in many cases, because if the data is nicely separable in these hyperplanes, then if you were going to use something like a random forest, you have to make cuts which are effectively always perpendicular or parallel to one of your axes. And so there, a random forest or even just a decision tree would have to make a cut that's here, and then on a subsequent run, it might have to take another cut in another dimension, just to build up this — what is a very easy math to figure out what these hyperplanes are. So this, again, is an optimization where you're trying to find all the hyperplanes that give you the maximal separation here, subject to some constraints: that if you had, let's say, a red value right here, the constraint would be, how many values of one class are allowed to be on essentially the wrong side of this line? So is that clear conceptually what's happening? There's some very nice properties of these things for these linear separator models, where you wind up getting some very, very good classification results.
Here's an example, though, where you've got data that might be hard to see, where it's blue here inside, and then it's red here in these two dimensions. And so you basically can't build a hyperplane that winds up separating these two things out. But instead, what you can do — if you click on this link here — oh, this link. Go to page. What you can do is actually project your data into a higher dimension than you had before, and then you can create hyperplanes that actually do give you perfect separability. So here, we're taking two dimensions. We're going to multiply those and push those into a third dimension. And now you can create a hyperplane that optimally splits between the blue and the red. And then you can project that back down into the dimensions that you had before. So it turns out that in a very, very large — in fact, infinite — dimensional space, you can build a perfect classifier, even for highly nonlinear relationships between the data.
The problem with that is you wind up having very, very large models. And again, you come back to this problem of the metric space that you're doing all of this work in. If this direction here is magnitudes and this direction here is which telescope you used, what is the metric in that space that you would actually wind up using? Again, I come back to things like random forest, because when you make your cuts, you're making it in the metric space of the thing you're cutting inside of. So it tends to be much more powerful.
All right. Let's go back to classification. Any questions about support vector machines? What I will say — it's not even anecdotally, I guess you could look it up — is that for lots and lots of different types of problems, both in astronomy and outside of astronomy, support vector machines rarely do better than things like random forest on real-world competition. They have other issues as well, like their model sizes can be very big. So if you try to move them around in memory and you have to use them often, you have to hand them to people, the typical rule of thumb is that they're about a third the size of the original data set, which can be quite large. Anyway, lot of people use support vector machines, but that tends to be mostly on the academic side. The other major issue with support vector machines is that there's no natural way to get probabilities out. Whereas with k-nearest neighbors or with random forests, you can just look at what all the different, effectively, votes are and say, well, I voted mostly for green, so that's the answer, but I had a couple other red classes for this new instance — then you can create a realistic notion of probability. With support vector machines, there is no way of getting a natural probability out of that. And many times, we want to get some notion of, I belong to this class with this probability. Yeah.
STUDENT: I think one of the main reasons that academics love support vector machines is that the optimization it does inside is closed form, in the same way that PCA is closed form, or the linear fitting, I think yesterday, is a closed form. And so it's truly a convex problem, closed form solution. So you know you're getting the best possible answer, even in the infinite dimensional reduction version. So it is a very beautiful algorithm, and it deserves its props there in that respect, the same way that PCA and linear fitting are just wonderful.
JOSH BLOOM: Yeah, again, I completely agree, obviously. And it's beautiful at the mathematical level. It just isn't practical. For all the reasons that we have up here, they actually can be very computationally expensive — not to build it necessarily, but to move it around. If you're doing a nonlinear kernel, it can be pretty expensive. But it's not all that informative either. You're making all these massive hyperparameter splits — it's not really clear what's going on. And then again, because we care about probabilities, support vector machines don't give you probabilities. And for me, the biggest knock is our data is heterogeneous in terms of what the different units are in the features, and we also have missing data. And so you can try to coerce what is a real-world data set or astronomy data set into something that you can build a support vector machine against. But in the end, you're basically doing a square hole in a round peg, whereas there are other models that allow you to natively take care of things like missing data.
And the last one is that in context of classification, support vector machines works in a 0 to 1 way. So if you're trying to build a multi-class classifier, you're effectively building up, do you belong to this class or another class? Do you belong to class B or the other classes? And then you wind up choosing the one that gave you the best answer, and that's your classifier. So it's not very natural for multi-class problems either.
STUDENT: I just want to add one more thing, which is that — this connection — all of the reasons that you don't like support vector machines are all of the reasons that I don't like PCA. And the two algorithms are very closely related, because they both are these closed form linear algebra solutions. And the reason that they work is because they depend on there being a well-defined metric and no missing data and homoskedastic and all that stuff. And so basically, whenever an algorithm is absolutely trivial, it's probably not appropriate for us.
JOSH BLOOM: For those that didn't hear that at home, this is a tweetable moment from Dave. If it's tractable — can I rephrase it a little bit? If it's tractable and it's beautiful mathematically, it's probably not that practical. OK. Now we're going to build a star-galaxy-quasar separator. Star-galaxy should be pretty easy if we had shape parameters. But again, we're going to use the original data we had, which is just photometry parameters. So this is a somewhat non-trivial problem. We're going to pull over 1,000 quasars, 1,000 stars, 1,000 galaxies. We're going to remove all the bad data. And we're going to take a look at our data. So here, we've got our object IDs and our different colors here. And we'll take a look at our different classes. So we've now got 3,000. This is what's called a balanced label problem — we've got one third, one third, one third. So this means the classifier really is going to have effectively an equal shot at making predictions on all of them. This obviously doesn't take into account that there are far more stars that you're going to be able to see in Sloan than you are going to be able to see quasars.
So again, we're not using our intrinsic knowledge about this. This is saying, given a source and given no other priors on what that is — if I told you I found a source in Sloan, the first thing you'd say is that's probably a star, because there are more stars in Sloan than other types of objects. If I said, well, it's not a star, you'd say, well, it's probably a galaxy. And say, OK, well, it's not a galaxy — well, then it's a quasar. So without giving you any other data other than saying it's not a star, you can pretty much figure out what something is. But given a new object in a purely frequentist way, without any priors, I want to know, given its data, what type of object is it? And here, we'll make our y value. So now, instead of it being redshift — a continuous variable — we're going to wind up having a three-class problem. And we'll run random forests on that with 200 trees.
And interestingly, before we see the results, you can actually count up the number of times that a feature is split upon during this process. And the number of times a feature is split upon is actually indicative of how important it is for the answer that you care about. Because I didn't go into all the details of exactly how you choose which feature to split on at every one of these different nodes — but effectively, there's a hyperparameter which says, choose a certain number of features to try out every time you get to a node. Figure out the one that is best at doing the split — this is what's called a greedy algorithm — and then split that way on that one. And the next time, just randomly choose another number. Typically, the rule of thumb is that you'll choose a number of features which is the square root of the number of features that you have. So if you have nine features, like we do, it will choose, at any split point, three random features to try. It'll figure out the optimal split across all three of them, and it'll choose the one that gives the best information gain. And then it'll split there on that feature. And then it keeps on going. So if you keep splitting on the same feature to improve your information at each of the nodes, it stands to reason that that would be more important. And there's a more rigorous way of keeping track of importance in your model with random forest. And then you can plot that up. Was there a question?
STUDENT: Yeah, how do you pick the number of estimators? What's the best–
JOSH BLOOM: So the number of estimators, this thing of how many features do I try — that's called mtry. How deep is the forest? How many instances should be in each node, in a terminal node? All of those things are called hyperparameters of random forest. And support vector machines is a different set of hyperparameters, et cetera, et cetera. We're going to show you how you do that sort of selection. But effectively, you have a bunch of rules of thumb of how many you should have, and then you wind up looking around that. Effectively, what you do is you either do a grid search or a random search over those hyperparameters. Again, with the little bird on your shoulder or devil on your shoulder saying, if you do this, you may be overfitting, you may be overfitting — so protect yourself. So there are ways to try to do that. But the answer is, logistically, you do a search over the hyperparameters.
OK, so u minus g color is pretty indicative of whether you're a star, a galaxy, or a quasar. The difference between your g-band magnitude and your aperture or your Petrosian magnitude is the next most important. These other two are important. And it turns out your i minus z color is not important at all, or not very important at all.
We can get what's called an out-of-bag error when we do our score, when we do random forest. Let's try to figure out what out-of-bag actually means. The way to think about out-of-bag is that whenever you are deciding to build up a tree, the first thing you do is you say, OK, I'm going to try a bunch of different features. But I also have to know all the different instances that are going to be part of my selection. You randomly choose your instances from the data, but you do it with replacement. So that means if you happen to randomly choose instance two, which could be, let's say, a star, you put it back into the bag. And then you choose again, and you might actually get number two. But you choose the total number of instances that you have. I'll do it in the column way. So if this is instances — it's the number of rows, and columns are the number of features — you're choosing, let's say, n instances. But you can do the math and figure out that you typically will randomly leave behind about 30% of your data at every single one of the top splits that you wind up doing when you build a new tree. And that's what's called the out-of-bag data. So the out-of-bag data, when you're building up a tree, is actually left behind. It is not part of building up the tree. But you can take out the things that you hadn't used, apply it, get what your answer is, and if you store up all the results, you get, effectively, the train-test split for you very nicely. So the out-of-bag score is, for all the data that randomly wasn't used, what were my predictions, and how well did I do? So this means that I got 95% of them correct with this classifier.
We can do the same thing in support vector machines, and there's lots of different kernels you can use to build up your SVM. And we'll print out what those look like and what the results are, hopefully. Any questions while this is running about the statements I just made about out-of-bag error? OK. So what am I doing here? Let's take a look. I'm just creating a bunch of different support vector machine models. And I'm creating different types of kernels, which is whether you're in a linear kernel or a non-linear kernel, which actually allows you to morph a little bit what the metric space is that you're using. And then there's some other hyperparameters of that model. Again, we're not going to go into the details of what all these different hyperparameters are, but it's useful for you to look at that if you want to play with it.
I don't know why this is taking so long. If this takes too long, I'm just going to jump back to my laptop. Just curious, what do we have here? Does anyone know how to figure out what kind of machine we're running on inside of a Unix? No? Is it cat etc proc-something? No. Cat /cpuinfo. Oh, it looks like a lot of cores. All right, let's see if this finished. What the hell? It's still going. All right, why don't I talk about hyperparameter optimization while this is running? If this takes another minute or two, I'll just jump back to my local version on my machine.
So this is where you would actually do this hyperparameter optimization that you were asking about. It really shouldn't take this long. Here, we grabbed this thing called grid search. And we can now create a dictionary that has a set of parameters, where we're now not individually creating all of those support vector machines by hand. We're going to create a set of parameters which will be fed in when we're creating our different models that we want to loop over. Or in this case, we're going to create a large-dimensional grid over all these three different parameters: what kernel it is, whatever gamma means, and whatever C means. I'm going to just kill this. Oh, we got some of the images. All right. So here's two of the support vector machines in the space that they wind up creating. Here, you notice that the linear kernel is splitting up the space as best as it can across the three different classes. Here, with radial basis function kernel, you wind up having some nice behavior here. Again, we don't know whether we're overfitting on this one blue dot or not. But that's just to show you what the different spaces are that get built out.
All right. So let's take a look at what we're doing here. GridSearchCV will run the fitter on this type of model, where it's passing different parameters and different combinations of those parameters. So we can run this here and now run it in a parallel way. Hopefully, this will be faster. And we can figure out what the best score is. So we're doing 168 different fits. That finished pretty quickly. So it looks like I have a lot of cores, but terrible CPU or something. And I got my best answer out here. And so you see, this actually got comparable to random forest. And here's the best model. The best model, in this case, was gamma of 0.1. It chose a radial basis function kernel or the other parameters we searched over. I think we looked at degree. So that's actually kind of nice, because we can run it over multiple cores. And there's something called joblib, which is what scikit-learn uses to actually push this out to as many cores as it can get access to. And it's running all of these embarrassingly parallel jobs for you. So we did 168 of those in seven seconds over however many cores there were here.
All right. So we can look at what's called the confusion matrix, which is a nice way of showing you how we wind up labeling correctly or mislabeling between classes 0, 1, and 2. I forgot what we called those things — 0, 1, and 2. So 0 is quasar, 1 is star, 2 is galaxy. And we can get an actual numpy array showing us what those values are. So this means that for all the things in the training set, we built a classifier, and then we applied it to the testing set. We got 231 of those right of quasar, whatever, 254 of star, and 255 of galaxy. And you see, these off-diagonals is where we had some confusion between those different classes.
If you get a perfect classifier and you see a confusion matrix that has no power off the diagonal and all the power in the diagonal, you've overfit. I guarantee it. Or you're using a trivial data set, where you shouldn't even be playing with it. But if you don't see some off-diagonal power, it means that something in your training set is leaking into your testing set, or something in your features knows about the answer in a way that you don't want it to know. If I had a label that said the label q, s, and whatever, g, probably we'd have a perfect classifier, because it would wind up figuring out that whatever column I wound up using in my feature data was actually exactly predictive of the thing that I wanted. Another way to figure out whether your classifier is any good is put the actual answers into your feature set. Run your classifier. And the most important feature better be the thing that's the answer. If it's not, you have a bad classifier. You've done something else horribly wrong. These are the types of things that you have to start thinking about as you're building these up.
OK, so let's do the same thing, where the question was how many trees do we use? We'll actually build one where we're looping over a whole bunch of different questions that we have for the hyperparameters of this random forest. So we're going to loop over 108 different types of models with different hyperparameters. Hopefully this will also finish in a finite amount of time. OK, that was pretty good. There we go. All right, so we got an answer, which is a little bit worse than what we had had before with the models that I chose before. We don't know whether this is an inferior one. It just happens randomly that we wind up getting slightly worse results, because again, there's a random component to the construction of this model. I thought I saw a hand somewhere. OK, no.
So here's our answer. Max_features 3, min_samples in the leaf node is 2. Number of estimators, 200. We can get those parameters, and we could save those for later. These are the best parameters of the random forest we built. Here's the best parameters of the support vector machines that we built. And I won't go into the details of what you can do in here, but you can actually set what your different scoring functions are. So if you want to do your correct choosing of which hyperparameters are best, you don't have to use whatever the default is, which is probably mean squared error for regression, and it's probably just accuracy or score for classification. Any questions about that?
It should be obvious that when you're doing grid search, it's not RAM-optimized and not compute-optimized. If you actually look at the code that's doing the grid search, it's doing a whole bunch of nested for loops. And oftentimes, it will wind up not saving pre-processing that you might do in the step. So if you had a pipeline that was part of your grid search, you might actually want to save the data and then reuse it for something that's inside of a for loop. These are not always optimal. When you have a lot of compute power and a lot of time to sit around, that's fine. What it turns out is that somebody wrote a paper fairly recently showing that if you just did a random search over your space, you get an answer which is pretty good, and it's pretty close to the right answer, because over multiple hyperparameters, it turns out most models are only sensitive to a few of those hyperparameters. And so searching over an entire grid in lots of the hyperparameters that you think you care about, you're spending and wasting a lot of computational time in a space that's not actually all that useful. There's also a whole Bayesian formalism for deciding what next hyperparameter to wind up using, given the results of all the other runs that you've had. That works irrespective of the model that you're actually building.
I want to jump back into the notebook and go into a little bit more details here, and some questions that you would be asking when you're doing hyperparameter optimization. How do I choose a model is really the critical one once you've done all the featurization and pre-processing steps. K-nearest neighbors is, what's the number of neighbors? Support vector machines is, what's the kernel that I'm going to use? What's the bandwidth? Random forest is, how many trees? What's mtry? With Gaussian processes, it's different sets of questions. I showed you this already. And then there's lots and lots of different metrics that you would potentially wind up using to create your optimization, and lots of different ways for you to show plots of how well you did on the models that you wind up choosing. Obviously, I won't go into all the details here. I will post this notebook in the Astro Hack Week GitHub repo, so you can see that for later on if you want to go into that. So all these things are available to you. And you may decide none of these metrics are the ones that I actually care about, because as we talked about later, we might not want to be optimizing on something like accuracy. We might want to be optimizing on something which is much more closely related to the problem that I have at hand.
What I will say — and this is from a paper we wrote in an astronomy context, where we looked at a bunch of image differences on the Palomar Transient Factory, and we had a bunch of labels of whether this thing that looked like it was a new source in the image difference was a real source or whether it was not a real source, as in a transient or not a transient. So this is what we call the real-bogus problem. And given the same input featurized data — across, I think it was 42 dimensions in this case — what you wind up seeing is that in all metrics for false positive and false negative, random forests wind up beating out everything else. And this is with the best tuning of each one of those different models. I won't claim that this is a universal type of curve that you'll see. But if you're asking the question, which model do I try, I often would try random forest first. And then maybe if that doesn't give you the results you want, you might do support vector machines. But the other ones tend not to do as well. And again, we already talked about this question about what we can pull probabilities out. So in this case, we want to minimize the false positive rate and minimize the false negative rate, the type 1 and type 2 errors. So you want to push that curve all the way back down to the bottom left. And at all places, we found that random forest actually did better.
We're not going to do the breakout session. I will continue to extol the virtues of random forest and say that random forest is built into Kinect. All this is doing inside of Kinect in hardware is running a random forest that's been pre-trained on lots of probably postdocs at Microsoft — maybe grad students, probably postdocs — where they walked around, and they then labeled which part of your body each pixel was. And then random forest is trying to figure out and colorize, in this case here, what body part you are in a given pixel. And so all it's doing is, at pixel level, it's saying what body part are you? What body part are you, given the input images? And random forest, because there's a bunch of if-then statements effectively, can run incredibly fast on the prediction side. So this is what they did. They took all this input data — left hand, right hand, shoulder, neck. They built a straw man model of what that would look like to create some notion of tracking. And then they just built a classifier at the pixel level. So here are the postdocs doing crazy things that you do in Kinect. They had a million training examples.
All right. So let me just finish on the hyperparameter optimization stuff. And then we'll jump back, and I'll show you some other types of models, show you some other examples of using machine learning in astronomy. All right. We're going to run this without doing a full grid search. And we got an answer which is pretty comparable to what I had before. We've been using this thing called joblib. Has anyone heard of Dask? OK, so a few people. Again, it's the people on the right. So it's the Python 3 people. Interesting. People at Continuum have built essentially a distributed computation and scheduler that allows you to make use of multiple threads and multiple cores and even multiple computers, where it keeps track of what the computations are that are going to be needed and are then going to be needed farther down in the computation tree. And it keeps in memory the things that it believes it's going to need later on. And then once it's done with that, it will wind up excising it out of memory. So it's a computationally and RAM-efficient way of doing computation. If you run install Dask and distributed, we can try to run this grid search using Dask, which will do the distributed thing in a completely different way. So here, I'm just going to pull over a very, very not safe for work — well, not safe for work in the traditional NSFW; not safe for work in the sense that it's not ready for production — way of this thing called dask-learn, which is trying to use Dask and rewrite all the scikit-learn fitters and optimizers. So we're just going to do that and pull that into our namespace. Hopefully that will work. Oh, nothing called Dask. I forgot to run this. Oh, right, because I'm running it in the cloud. I'm not running it on my machine. Installing Dask. Extracting packages, complete. OK, so this should now work. Yep. That worked, right? Yeah.
So now, we'll do the same exact thing. But instead of doing the grid search that comes with scikit-learn, we'll use a different grid search. And this one actually may take a little bit longer. We'll see how long it takes, just because it sets up a lot more infrastructure. But if we did many, many more parameters, or we did this on a very large cluster — let's say it took five minutes to build one of these models — this would totally win. So we'll let that run. I won't execute this thing here, which allows you to build up a cluster and then run it over a cluster. I'll come back to a different kind of clustering. So here, this took 31 seconds. We got about the same answer as we had before. Let's see how it did timing-wise relative to the other one. That one took nine seconds. So still much faster on the scikit-learn side. But there are workloads that you can find on the web using dask-learn that actually do much faster and much more computationally efficient. What I didn't do is keep track of the RAM usage. I think one of the conceits of Dask is that it's just going to be more efficient at using RAM. All right. I'm going to jump back now to a couple other things I wanted to show you. Who's heard of deep learning? OK, now almost everyone, not just Python 3 people. Effectively, the idea around deep learning, or what people are very excited about, is not only does it give very, very good and probably the most accurate answers on a class of problems — typically around language understanding and around voice-to-text, around image processing, around video processing — the thing that people really like about it is that you don't have to do a lot of featurization. In fact, the idea is you can throw raw data at a deep learning network, which, in the end, is just taking your raw data, and it's deciding how to combine different pixels from different places on, let's say, an image, with weights that it winds up learning through a bunch of very clever techniques that allow it to learn in a finite amount of time. And then it takes the results of whatever it wound up taking from your raw data. That's going into what's called the next layer. And then you have another set of weights which are multiplied against that data to get to another layer.
And the idea for image processing is that you don't have to know a lot about what it means to be an edge to detect something that is an edge. You don't have to know what it means to be a cat to find cats in images. As long as you have an objective function on the other side that you're trying to, let's say, classify am I a person or a cat — if you're trying to get a two-class problem as an output — if you have enough data, you let the machine figure that out. So what's happening — unlike what we just showed you throughout this whole tutorial, where we took our raw data, we did some feature engineering, we threw some stuff out, we built some new features, and then we did our classifier — in this case, with deep learning, the idea would be that you can really do both the feature engineering and the classification in the same sort of process.
This is very computationally expensive, although again, there have been a large number of advances in the field that make it at least tractable. But oftentimes, because a lot of this is a lot of matrix multiplication, these are things that you can actually wind up doing on GPUs very effectively, or even FPGAs. It's worth taking a look at this result here, where people built up a network to do some inference on galaxy images, essentially to classify galaxy images, using raw Sloan data. And here, you can see that their input was RGB data. And they wind up building up a network to get to some sort of answer that they wind up caring about. I'm not going to go into the details of what all these different types of layers are and how you wind up combining all of the data, other than to say that yes, you give up on feature engineering. That's a big plus one for deep learning. But in the end, there's this whole black art of how you build up your network. And how did this person get to the point of building a network that looks like this? Because everyone's network is going to wind up being different. And you wind up getting different results depending upon what your network is and all the different hyperparameters of the training of that network. All of that is stuff you still have to do.
So we've pushed featurization, which, as domain experts, I tend to like a lot, because it allows you to commune with your data. It allows you to ask questions that you know are physically sound and relevant for the question that you're asking of your data. And it allows you to build features that are going to be informative and then throw it into a classifier. Whereas in this case, you wind up potentially just throwing everything in, the whole kitchen sink, and hope that you get a good answer out. Now, there's nothing forbidding you from using deep learning networks to work on featurized data. But because these tend to work well in the context of two-dimensional data, most of the time, you're not doing featurization. And then there's some sort of notion of some local connection between pixels. In that case, you certainly could build a deep learning network on an input vector, which is what we've been operating on now. But it typically will wind up wanting to work on two-dimensional images. Or at least that's how it's been used to most effect recently. Any questions about deep learning? Yes.
STUDENT: I would just say that that blog post by (INAUDIBLE) here is excellent, and it's really worth reading. And if you're interested in thinking about deep learning and the surrounding context, it covers a lot of the ground. (INAUDIBLE) choices in there, but it's a very good, straightforward introduction to how you (INAUDIBLE).
JOSH BLOOM: Just like machine learning is not the answer to all your inference problems in astronomy, deep learning is not your answer to all machine learning problems in astronomy. And in fact, in the image problem, where we looked at the real-bogus, we did a traditional featurization step on that, and then we did random forest. We also tried using deep learning, and we got inferior answers. And the reason being is that while you can take networks that have been pre-trained on other images and apply it to your image, what's actually happening in many of these cases is it's learning a lot about the details of the questions that you're asking of the original image set. And we didn't have enough data to train what is effectively millions and millions of nodes in this network. And so we just got inferior answers. So here's an example where we needed to put it in production to work on real-world telescope data effectively in real time, and deep learning just didn't work. You could argue we're not deep learning experts and we didn't know how to build the networks correctly, but we did spend enough time on it. I think Danny kind of banged his head on that enough. But it's not going to solve all of your problems.
Another approach of using deep learning, which I think is gaining more and more traction, is to use this notion of autoencoding, for those that have heard of that, where you don't try to predict an outcome with labels. You try to build up the data that you had before. So you have a deep learning network that looks like this. And what you're doing is you're compactifying your data and getting it more and more summarized over time, so that by the time you get to these deep layers, you, in principle, have learned concepts about your data. Whereas at the very high layers, especially for image processing, you're, in principle, learning sort of low-level features that you would apply to data if you were going to do it from scratch. So this could be like histogram of gradients. It could be other types of low-level featurization that you would ordinarily do if you didn't know about deep learning. It's been shown that some of these higher-level layers closer to the raw data actually do some traditional filter bank techniques and are learning those intrinsically, which is pretty amazing. But farther down, you've got concepts and deeper concepts of what it means to be the data that you've given it.
If you then take this compact notion of what it means to be your data, and then you actually start building and start expanding that data back out, you can try to build something that doesn't get smaller and smaller, but now gets bigger and bigger. And then you wind up trying to predict the data that you started off with. That's what's called autoencoding. What's nice about autoencoding is that you don't have to come at this problem with the classification set in mind. You can just say, I've got a whole bunch of data — a bunch of images, let's say. I want to build a deep learning network that gives me back the data I started off with. But now, I can use this stuff in the middle, which are the deep concepts, and I can, in principle, use those as features for another learner. So what we've been doing in my company, what we've started thinking about in the astronomy context, is marrying both deep learning and something like random forest, where the deep learner is not the classifier itself, but is just the thing that's building the features for us. That's definitely something worth exploring. If somebody wanted to build a hack around that this afternoon, that'd be pretty exciting. Happy to talk with them about it.
The other thing I should say, in the context of deep learning and more generally about machine learning, is that when you have more data, you wind up traditionally doing better. So even if you have the same learner, you have the same featurization process. And this comes — something Peter Norvig at Google said a while back: more data beats clever algorithms, but better data beats more data. So get more data, but make sure it's pretty good. And oftentimes, if you want to improve your model — given that you've just spent a bunch of time featurizing it and extracting what you think is the most information out of this and then throwing it into a good classifier — the only way to get a much better result is typically getting more data. There was a question.
STUDENT: Yeah, the question was related to this. How good example size is good enough for deep learning? How many tuples is good enough to move to a deep learning model? Especially, is there a rule of thumb between the number of features that you have and the number of examples?
JOSH BLOOM: Oh, good. So the question is, what's the rule of thumb around whether you do deep learning or something else? I don't have a good one. What I will say is a couple things. One is, if they're images or they're sonograms of frequencies or something like that, that may be a pretty interesting example of where you'd start thinking about using deep learning, because those are the obvious places where they're doing really well. If you have lots of metadata and heterogeneous data, in the sense that you have categories of your features, and you also have numerical values and maybe even strings, that's a pretty terrible place to think about deep learning. And the only time you'd have to do it is after you did a whole bunch of featurization. But once you've done all the featurization, you might as well just use a traditional learner. So if your data looks like that, that's a little scary for deep learning. If you're looking at pixel data, and every pixel goes from 0 to 256, that's a good place to think about it. Now, in terms of the volume of the data, millions of images. You can do well with, say, tens of thousands of images. There's a whole notion of what's called transfer learning, which I alluded to, where you could take a network that's been built on images taken off the web, and then take that whole network and just use it to apply it to your images. And in principle, there are features that it's effectively learned and concepts it's learned that it could wind up using and reusing, that could be informative and allow you to use smaller amounts of data. I think once you're into the millions of instance level, then you can credibly start thinking about deep learning. People on the web may be flaming me right now, because they say, no, there's an example where we can do 10. If you're at the 10-instance level or 100 or 1,000, you have to be using some of these other techniques. Any other questions?
STUDENT: You mentioned the concept of rights. (INAUDIBLE) Is this something that you have to set up (INAUDIBLE) yourself and then validate over there?
JOSH BLOOM: I missed the beginning of your question. I mentioned the concept of what? Rights?
STUDENT: Weights. Weights. (INAUDIBLE) the raw data and then the weights. You have the layers that you want to choose. (INAUDIBLE) I just wonder, when you said the machine or the algorithm learns later on what (INAUDIBLE) to choose, I wonder whether that's something that you said initially and then–
JOSH BLOOM: Right, so the process by building up these networks. What should be clear, just to step back, is that a big part of the work that people do in machine learning at the academic level — and a lot of that's moved into R&D centers inside of companies — is to figure out optimal ways to do the optimization problem of, how do I build up the network? Or in random forests, what's the optimal way to make the splits to get the answer out that I want, where optimal now can be not just accuracy, but RAM efficiency and CPU cost, et cetera, et cetera? So a lot of the work is on building up the networks and learning how to learn. Once the network's been built, then it's just a matter of just cranking data through. And one of the nice things about things like random forests and deep learning networks is that it's not that computationally expensive to take the data, throw it through, and you get your answer out really, really quickly with not a whole lot of extra math. But so the question about how the weights are created — that comes down to how you decide to build the network. Oftentimes, the weights will be randomly assigned at the beginning. And then you have this notion of what's called backpropagation, where you start from the answer. And there's a way, in a mathematically tractable way, to start and update the values of the weights going backwards, forward. There's also forward ways of doing it as well. So that's all of the techniques that people are working on to try to make that better.
All right. We already talked about improving models. I guess what I'll say — I'll just show you a couple examples of — how are we doing on time?
STUDENT: It's 11:36.
JOSH BLOOM: OK. Couple examples of machine learning approach to classification of what we've done. I told you about the 0, 1 classifier of real-bogus. We've also done this on variable stars. So now we've worked in the time domain, and we've taken photometry and other things, and we wind up building a classifier off of that. So here, we took what is very heterogeneous data. Sometimes we've got photometry in r-band magnitude. Sometimes it's regularly sampled data. Oftentimes, it's not regularly sampled data in time. And instead of throwing the raw data in, in RA and dec and all that stuff, we build a whole bunch of features off of that. So you can build features in the frequency domain. You can build features with unordered statistics, et cetera — so variability metrics, periodic metrics, shape metrics on the light curves, and then context metrics of where this thing is in the sky and what's near it. And we were able to get some pretty interesting results out of that.
We wound up building up — and this is now still a work in progress with some people that are in the room — essentially a framework that allows any of you all, and hopefully a lot of your colleagues, to start getting not just the pipeline that I showed you at the IPython notebook level, but even data handling and project handling and reproducibility notions out of your data. This is a project we call Cesium, and we just released the initial version of this. This provides a lot of the featurization capabilities around time series data, and then also a whole bunch of access to scikit-learn modeling and feature selection, et cetera, and then even plotting that allows you to build up entire frameworks even around some inference problems you might have in the context of time domain. This is different than scikit-learn in the sense that scikit-learn is a set of models, and it focuses very much on the modeling part, not so much on the featurization part. We've spent a lot of time on the engineering components around the featurization. And our idea is that we're going to be able to take this to other domains as well, like in seismology and neuroscience, and take what is fundamentally a similar type of signal — just measurements as a function of time to predict outcomes — and use essentially our feature bank to actually get good answers. So I'm happy to talk with people about that if they want to start using it and playing with it. There's a way to install it and get it working on your laptop with Docker. So we applied our classifier that we built on variable stars to something like 50,000 stars from ASAS, where we only had a training set of 810 over 25 different classes. So this was a hard problem — our training set size was really small. And we try to get, out of light curves like this, what the probability of it being an RR Lyrae is. And we try to do that as far back to the raw data as we could. And we got something like a 15% error rate, which is pretty good across these multiple different classes.
And we built a website called bigmacc.info, if people want to check it out, where it allows you to peruse through the hierarchy of variable stars and then click into those. So here's pulsating, and we'll go into different types of RR Lyrae. We'll get a fundamental overtone RR Lyrae out of this. And these are the ones that were predicted, 405 of them out of the ASAS catalog. And we show the probability of it belonging to that class. And then you see the probability vector of what comes out of that, of belonging to the class. You see the raw data and the folded data. And then we made it social so that Facebook could buy us one day.
But anyway, for those of you that are working in variable stars, the whole idea is that you've got to start getting used to probabilistic catalogs. So when we say this belongs to this class or that class or that class, the real answer is it belongs to this class with this probability, and it belongs to that class with that probability. And oftentimes, as astronomers, we want to take spectra of things that are of that type. But intrinsically, we know that there's some chance that they may not be part of that type. We're trying to formalize that in studies like this that allow us to put these different objects into different buckets and allow us to do science across those different buckets. I can go into the details of that if people want to ask me about it offline.
The last thing I'll say — and then I think I'll probably end, and leave for you in the notebook, you'll see there's some stuff on doing anomaly detection — is that we did a study with a student of mine named Adam Miller, where we looked at Stripe 82 in Sloan, which was effectively a five-color study in time of all the stars that had essentially shown some level of variability over several years. And then we built a classifier on the variability metrics, using all the parameters that I showed you before, to predict stellar quantities that you could only get out of spectra. So that's temperature, gravity — log g — and metallicity. And we got about 5,000 spectra or so ourselves, or they existed in the catalog already. And we wound up showing that using time domain parameters plus colors, we're able to get root mean squared errors which are comparable to what you would get out of low-resolution spectroscopy itself. So here, we used random forests. We did lots and lots of featurization, tried to protect ourselves as much as we could against overfitting. But this is just another example of using machine learning in real-world contexts. What we wouldn't claim is that you could take our model and then apply it to another time domain survey and get as good of an answer out of it. That's a whole separate problem. But in principle, if we got more Sloan data in other parts of the sky without taking spectra, we could, in principle, learn what the fundamental parameters are of that star or stellar system.
All right. I'm going to end there. We'll have time for a couple of questions. Let me just say, again, just to reiterate: think of machine learning as another set of tools in your toolbox. And if you haven't been trained on it — and by the way, I've just given you an opening to all of this if you haven't seen this before, so don't consider yourself well-trained — there are lots of ways to hammer your thumb with your new tool. And you may think you're building this amazing house, but in the end, you've built complete crap. That happens a lot. Machine learning is fraught with places where you're introducing biases that you didn't know about. But it can be very powerful. And as long as you're asking a question that's appropriate of the data that you either have or you plan to get, where machine learning is the right type of tool, you're off to a good start.
And when in doubt, it's worth asking other people in your field. And in particular, for those with home institutions that have stats departments and people working on machine learning, it's not a bad idea to go to them and say, here's the type of problem I'm trying to do. Give them a little tutorial for five minutes on the kind of science you're trying to accomplish. And let them tell you whether they think this is appropriate for machine learning or whether it's not. Because oftentimes, it's not appropriate. So just bear that in mind.
I often come back to this nice quote that Jim Gray — who was at Berkeley and then at Microsoft, sort of prototype of a modern data scientist — whose quote was, “I love working with astronomers because their data is useless.” And he meant it, because if you're in Microsoft and you're trying to build new algorithms against data, the data that you typically have access to has personally identifiable information. If it leaks out that you're using it, it's really bad. But he loved working with astronomers because our data — who cares? Like, oh wow, you showed me a galaxy that I wasn't supposed to see. It's not the end of the world. You're not going to start a war. So statisticians and computer scientists actually really like working with astronomers, because our data is pretty big. It has some interesting properties — noise aside; most of them don't like to think about that. But there's lots of questions to ask of that data, where machine learning may actually be useful. And they can hone models, try new scaling curves, get new benchmark data sets around that. But the flip side isn't always true. Just because you've got data and you're an astronomer and you've been at this Hack Day doesn't mean that you should be trying all of these out against all these cool, fancy new approaches that you've been reading about at blog post level. So just be careful. And if you're careful, I think you'll go a long way. So I'll end there, and happy to take a few questions. [APPLAUSE] Yes.
STUDENT: On that note, have you seen machine learning being used in contexts where you could tell straight away, that's completely inappropriate?
JOSH BLOOM: So the question is, have I seen machine learning being used in places where it's completely inappropriate? Yes. Oftentimes, if it's truly inappropriate, it won't make it all the way out to the public sphere. But there are lots of examples from people who reside in computer science and stats who have applied machine learning to questions that aren't all that interesting to astronomers, using astronomy data. I won't name names. But they have the other problem of not knowing the right questions to ask, and not being able to evaluate independently, other than at the accuracy level, how well they're doing. So if you build a photometric redshift estimator on photometry, and you only look at the r-band magnitude or something — even if you get a classifier which is pretty good, that's OK, but we know that for all intents and purposes, when you're doing photo-z's, you're doing it so that you can do something in cosmology. And there, you need to also have your uncertainties, but you probably also bring in lots of other data to bear. And so looking at this one benchmark data set that has just that data available just is cool because they use astronomy, but it's not all that useful for any of us. Any other questions? Dave.
STUDENT: Just wanted to make a comment about probabilistic catalogs. Maybe you thought about it a little bit. We're obviously very interested, in my group, in probabilistic catalogs. And one issue is that if you want to combine data from different sources, you really want likelihood information to combine data, not posterior information. Implicitly, random forest generates posterior probabilities that things are in different classes, and so there's an implicit prior that's been applied. And — just first of all, just as a general comment — it's very important that when we release probabilistic catalogs, we also release the prior information, the implicit prior information. But can you say a few words about what the implicit priors are for the random forest classifier?
JOSH BLOOM: So one obvious prior is that the distribution of your labels in the data you haven't applied your model to yet is the same distribution as in your model. If you're building a star-galaxy classifier, again, there's vastly more stars in the Sloan catalog than there are going to be quasars or galaxies — I think that's true. And so if you now apply it to just a pure, let's say, subset of galaxies, you're going to wind up getting answers that don't make sense, because I already knew that they were all galaxies. Why am I sometimes getting the answer of star? So that's there. And the imbalance in the data set, in principle, is being learned by all these different classifiers. You can figure out ways to, post-hoc, pull that assumption out without breaking the classifier itself. But that's an obvious one.
The other one, which is obvious, is that you're assuming that the classifier that you then apply on new data isn't just looking at the same universe. It's acquiring the data with the same noise properties, done in the same part of the sky. So again, if I built a star-galaxy classifier and I just happen to use it around places where there were a bunch of stars, like in the galactic plane, and then I applied that blindly to a place somewhere else in the sky where there's just more galaxies around, you're going to get wrong answers. Or if I did it not on the Sloan catalog, but now I did it on PTF, I'd get wrong answers. So there's the obvious mismatches of what you're actually looking at and how you obtain the data. There's probably other, more subtle priors that I'm not coming up with at the top of my head, but those are the obvious ones. OK, thanks.
STUDENT: Thanks again, Josh. [APPLAUSE]